PBR（基於物理的渲染）學習筆記2

木的樹發表於2021-04-22

理論

PBR是基於物理的渲染，核心是從能量守恆角度將各個方向的光源進行積分。BRDF就是根據物體的各種性質經過一系列實驗得到的一個雙向反射分佈函式來進行模擬。在之前那篇文章中，首先對於環境中已有的點光源和方向光源分別進行caculatorFinalColor處理根據diffuse和specular；

vec3 calculateFinalColor(PBRInfo pbrInputs, vec3 lightColor) {
  // Calculate the shading terms for the microfacet specular shading model
  vec3 F = specularReflection(pbrInputs);
  float G = geometricOcclusion(pbrInputs);
  float D = microfacetDistribution(pbrInputs);

  // Calculation of analytical lighting contribution
  vec3 diffuseContrib = (1.0 - F) * diffuse(pbrInputs);
  vec3 specContrib = F * G * D / (4.0 * pbrInputs.NdotL * pbrInputs.NdotV);
  // Obtain final intensity as reflectance (BRDF) scaled by the energy of the light (cosine law)
  return pbrInputs.NdotL * lightColor * (diffuseContrib + specContrib);
}

對於環境中光源的處理已經完成了，但是前面說到PBR是對各個方向的光進行積分，即對環境中各個方向能夠反射進入人眼中光都需要處理。

首先光怎麼來，可以認為是從環境貼圖中來，我們就認為環境貼圖中每個畫素顏色代表代表一個微分光源，也就是說要對環境貼圖中所有紋理進行遍歷求和，這個過程顯然對於實時渲染時不可接受的，這麼這個時候就出現了蒙特卡洛積分。蒙特卡洛積分的思想是在整個積分割槽間內，隨機的進行有限個取樣，通過取樣點的均值來進行近似。

想具體瞭解下蒙特卡洛積分的，可以看這篇文章，這是迄今我見過最通俗易懂的文章：https://blog.csdn.net/i_dovelemon/article/details/76286192

但是呢光有基礎蒙特卡羅並不行，因為我們是隨機取樣，每個取樣點實際上對於整體的貢獻度是不相同的，所以我們還需要計算每個取樣點對於整體的權重情況，那麼這個計算重要性的過程與蒙特卡洛結合就稱為基於重要性的蒙特卡洛積分https://blog.csdn.net/i_dovelemon/article/details/76786741。

現在整個公式變成這個樣子了

再來現在取樣有了，權重有了，還有一個問題要解決就是取樣的分佈情況。如果都集中在高權重或低權重對整體結果的影響是很大的，所以圖形學在這個過程中有專門的一個取樣序列的問題。通過一個有特點的取樣函式來生成一些取樣點，這裡使用的是Hammersley取樣序列演算法https://blog.csdn.net/i_dovelemon/article/details/76599923。

好了，下面來看下我們彙總後的程式碼（我們對於環境光來說漫反射部分實際上各個方向都是一樣的，所以這裡重點看鏡面反射部分）：

float3 SpecularIBL( float3 SpecularColor , float Roughness, float3 N, float3 V )
{
float3 SpecularLighting = 0;
const uint NumSamples = 1024; // 使用了1024個取樣點
for( uint i = 0; i < NumSamples; i++ )
{
float2 Xi = Hammersley( i, NumSamples ); // 計算一個隨機取樣序列

float3 H = ImportanceSampleGGX( Xi, Roughness, N ); // 將一個二維取樣序列轉換成三維空間中的取樣方向
 // 下面是計演算法線、取樣方向、視線等各種方向的一堆夾角
 // 看上圖L是從一個隨機取樣方向計算出得到的環境入射光源的反方向
float3 L = 2 * dot( V, H ) * H - V;
float NoV = saturate( dot( N, V ) );
float NoL = saturate( dot( N, L ) );
float NoH = saturate( dot( N, H ) );
float VoH = saturate( dot( V, H ) );
if( NoL > 0 )
{
    // 計算環境光源顏色，envMap很可能是立方體貼圖
float3 SampleColor = EnvMap.SampleLevel( EnvMapSampler , L, 0 ).rgb;
// 下面是計算BRDF的specular部分的G和F，這裡並沒有計算D，因為在BRDF/pdf過程中，D被消除掉了。
float G = G_Smith( Roughness, NoV, NoL );
float Fc = pow( 1 - VoH, 5 );
float3 F = (1 - Fc) * SpecularColor + Fc;
// Incident light = SampleColor * NoL
// Microfacet specular = D*G*F / (4*NoL*NoV)
// pdf = D * NoH / (4 * VoH)
// 上面pdf公式可以看出重要性跟粗糙度、法線與取樣方向、視線與取樣方向相關的。粗糙度是一個經過物理實驗測量的值
SpecularLighting += SampleColor * F * G * VoH / (NoH * NoV);
}
}
return SpecularLighting / NumSamples; // 求和再取均值就是蒙特卡羅積分的體現
}

上文中講到了先有采樣序列，然後是將一個二維隨機數對映為取樣方向，下面就來看下這個過程：

要理解這個，得要有立體角的概念，這個概念如果不明白可以搜一搜，講的挺多的。

對於一個微分立體角來說，要確定一個向量只需要有兩個量，phi和theta；這兩個量就是通過上文中的取樣序列生成的。下面要做的就是把這微分立體角座標轉換成三維空間座標。（注意：圖片中說的鏡面反射方向應該是法向量方向，這裡可能是GPU GEM中文作者翻譯錯誤了）

float3 ImportanceSampleGGX( float2 Xi, float Roughness, float3 N )
{
float a = Roughness * Roughness;
float Phi = 2 * PI * Xi.x; // 水平方向的phi
// theta不知道是怎麼計算出來的，可能也是根據一個數學理論來計算的，這裡可以看到有將粗糙度考慮進去
float CosTheta = sqrt( (1 - Xi.y) / ( 1 + (a*a - 1) * Xi.y ) );
float SinTheta = sqrt( 1 - CosTheta * CosTheta );
float3 H;
// 根據微分立體角座標求得以該表面為原點，鏡面反射方向為微分球的區域性三維座標系
H.x = SinTheta * cos( Phi );
H.y = SinTheta * sin( Phi );
H.z = CosTheta;
// 求表面切空間的的基底，切空間基底向量座標系為世界座標系
float3 UpVector = abs(N.z) < 0.999 ? float3(0,0,1) : float3(1,0,0);
float3 TangentX = normalize( cross( UpVector, N ) );
float3 TangentY = cross( N, TangentX );
// Tangent to world space
// 下面這個操作的前提是需要微分球的縱軸方向要與該點法向量的軸重合才行，而不是圖片中說的鏡面反射方向，這裡可能是GPU GEM中文作者翻譯錯誤了
// 因為微分球轉換的三維座標與且空間重合，所以這裡是將微分三維座標進行向量分解，最終得到一個三維空間座標下的單位向量
return TangentX * H.x + TangentY * H.y + N * H.z;
}

可以看到這是在實時渲染情況下的計算過程，根據GPU 精粹3中第20章的結果來看，並不用1024次取樣只需要40幾次即可。

這種實時渲染的優點真的實時渲染，不需要提前生成BRDF 的LUT查詢表，但問題也是每幀都這麼計算還是很耗費效能，所以後來的大牛各種研究，就是我們在PBR學習筆記中看到的那部分，在IBL這部分直接讀取紋理。

其實在圖形學繼續深入的研究方向比如模擬、光線追蹤，大部分處理過程都是一個近似過程，這也是為什麼像GPU精粹、GPU Pro、GPU Zen中經常開頭一複雜積分，到了最後程式碼過程其實還是簡單的加減乘除和迴圈。所以在整個模擬的過程中，首先通過物理理論研究完成完整計算公式，後續進行一步步近似和簡化，逐漸轉化成GPU可執行的程式碼程式。這也是圖形學奇高的門檻，一個數學公式讓人望而生畏，轉頭就跑。好了下面說一下怎麼進行進一步的近似處理。

經過一系列數學研究之後拆成了兩個公式。拆成這兩部分的優勢是，兩部分都可以做預處理，通過程式提前寫入紋理中，實時渲染只需要去紋理中去查詢即可，做簡單的加減乘除。

前一部分可以預處理成跟roughness和cos(v)相關的光照顏色，原理跟上面實時渲染差不多，取了取樣了一批環境紋理的顏色，取均值。

float3 PrefilterEnvMap( float Roughness, float3 R )
{
float3 N = R;
float3 V = R;
float3 PrefilteredColor = 0;
const uint NumSamples = 1024;
for( uint i = 0; i < NumSamples; i++ )
{
float2 Xi = Hammersley( i, NumSamples );
float3 H = ImportanceSampleGGX( Xi, Roughness, N );
float3 L = 2 * dot( V, H ) * H - V;
float NoL = saturate( dot( N, L ) );
if( NoL > 0 )
{
PrefilteredColor += EnvMap.SampleLevel( EnvMapSampler , L, 0 ).rgb * NoL;
TotalWeight += NoL;
}
}
return PrefilteredColor / TotalWeight;
}

第二部分，最終而後半部分可以轉化為一個kx+b的線性函式，這裡，只有x未知（可以認為這裡的x是物體的specularColor），而k，b只跟roughness和cos(v)有關。而且roughness和cos(v)都在[0,1]之間，所以，我們可以生產一張如下的紋理，建立一張對映表，通過roughness和cos直接找到對應的k，b值，省去中間大量的取樣計算。當然，紋理越大，越精確，但畢竟只是[0,1]之間的插值，所以是近似值。

雖然大部分地方顯示的是上面的圖，但是在這個庫中https://github.com/KhronosGroup/glTF-Sample-Viewer/tree/glTF-WebGL-PBR，使用的是下面的圖，看到y軸剛好反了下：https://github.com/KhronosGroup/glTF-Sample-Viewer/blob/master/assets/images/lut_ggx.png

而在上一期文章中，獲取brdf的程式碼中給紋理的座標也確實做了1減的處理

好，回過頭來，根據UE4的方案程式碼如下：

float2 IntegrateBRDF( float Roughness, float NoV )
{
float3 V;
V.x = sqrt( 1.0f - NoV * NoV ); // sin
V.y = 0;
V.z = NoV; // cos
float A = 0;
float B = 0;
const uint NumSamples = 1024;
for( uint i = 0; i < NumSamples; i++ )
{
float2 Xi = Hammersley( i, NumSamples );
float3 H = ImportanceSampleGGX( Xi, Roughness, N );
float3 L = 2 * dot( V, H ) * H - V;
float NoL = saturate( L.z );
float NoH = saturate( H.z );
float VoH = saturate( dot( V, H ) );
if( NoL > 0 )
{
float G = G_Smith( Roughness, NoV, NoL );
float G_Vis = G * VoH / (NoH * NoV);
float Fc = pow( 1 - VoH, 5 );
A += (1 - Fc) * G_Vis;
B += Fc * G_Vis;
}
}
return float2( A, B ) / NumSamples;
}

這裡有個庫，可以用來生成這種查詢表：https://github.com/HectorMF/BRDFGenerator

那麼按照UE4的論文，根據近似公式之後，這部分IBL的著色器程式碼變為：

float3 ApproximateSpecularIBL( float3 SpecularColor , float Roughness, float3 N, float3 V )
{
float NoV = saturate( dot( N, V ) );
float3 R = 2 * dot( V, N ) * N - V;
float3 PrefilteredColor = PrefilterEnvMap( Roughness, R );
float2 EnvBRDF = IntegrateBRDF( Roughness, NoV );
return PrefilteredColor * ( SpecularColor * EnvBRDF.x + EnvBRDF.y );
}

可以看到完美對應近似公式，上面說的第二部分拆成了kx+b的形式