从CubeMap生成Lambert下的环境辐射率图

前言

本篇主要是对《Coding Labs :: Physically Based Rendering》中内容的总结与实践。

本篇将对一些原理进行概况,但可能理解并不准确。
随后,将实践原文中的代码并观察效果。

渲染方程

概括上讲,渲染方程所描述的是:在知道一位置的所有入射光的情况下,这个位置某一方向的出射光是什么?
L o ( p , ω o ) = ∫ Ω f r ( p , ω i , ω o ) L i ( p , ω i ) ( n ⋅ ω i ) d ω i L_o(p,\omega_o) = \int_\Omega f_r(p,\omega_i,\omega_o)L_i(p,\omega_i)(\bf{n} \cdot \omega_i) \it{d} \omega_i Lo​(p,ωo​)=∫Ω​fr​(p,ωi​,ωo​)Li​(p,ωi​)(n⋅ωi​)dωi​

方程中基础符号的意思如下:

符号 意思
p p p 位置
n \bf{n} n 表面的法线方向(宏观)
i i i 与 o o o 入射 与 出射。
ω i \omega_i ωi​ 与 ω o \omega_o ωo​ 入射方向 与 出射方向。
L i L_i Li​ 与 L o L_o Lo​ 入射光 与 出射光。
f r f_r fr​ 双向反射分布函数(BRDF),他表示了一条入射光在各个方向上的出射情况。
∫ Ω \int_\Omega ∫Ω​ 在半球面上积分

那么,组合起来的意思就是:

符号 意思
L o ( p , ω o ) L_o(p,\omega_o) Lo​(p,ωo​) 在 p p p 位置上, ω o \omega_o ωo​ 方向上的出射光。可以看作是需要求出的结果。
L i ( p , ω i ) L_i(p,\omega_i) Li​(p,ωi​) 在 p p p 位置上, ω i \omega_i ωi​ 方向上的入射光。可以看作是需要提供的数据。
f r ( p , ω i , ω o ) f_r(p,\omega_i,\omega_o) fr​(p,ωi​,ωo​) 在 p p p 位置上, ω i \omega_i ωi​ 方向上的入射光有多少比例作用到 ω o \omega_o ωo​ 方向上。例如对于一个绝对的镜面反射比如镜子,它的情况会是这样:当 ω o \omega_o ωo​ 精确为 ω i \omega_i ωi​ 的反射方向时,函数结果为 1.0,否则函数结果就是 0.0 。
( n ⋅ ω i ) (\bf{n} \cdot \omega_i) (n⋅ωi​) 表面法线方向点乘入射光方向,也就是乘算上入射角的余弦。
∫ Ω f r ( p , ω i , ω o ) L i ( p , ω i ) ( n ⋅ ω i ) d ω i \int_\Omega f_r(p,\omega_i,\omega_o)L_i(p,\omega_i)(\bf{n} \cdot \omega_i) \it{d} \omega_i ∫Ω​fr​(p,ωi​,ωo​)Li​(p,ωi​)(n⋅ωi​)dωi​ 在整个半球上对所有入射光对 ω o \omega_o ωo​方向上的贡献进行积分。(为什么不是一个完整的球?因为不需要算背面半球的光)

如果这一函数能够解答,那么就可以求出场景中某一点到眼睛的方向的出射光,也就是应渲染出的颜色。

CubeMap

要计算渲染方程所遇到的第一个问题是:如何知道场景中一点的所有方向上的入射光?

CubeMap 是一种方案,其将所有方向的入射光存到六个面的立方体中。

在之前的博客《图形API学习工程(21):使用CubeMap纹理》中已经实践了相关操作。

Lambert

接下来遇到的问题是关于BRDF。

首先,BRDF肯定是取决与材质的。比如对于一个绝对的镜面反射,他的BRDF就很简单:只有在精确是反射的方向上BRDF结果为1,其余都是0。那么在这种情况下就只需要知道那一条完全反射的光线就可以了。

但是对于更粗糙的表面呢?四面八方的入射光都需要考虑,那么就需要进行实时的积分了,而这样的计算在实时中计算速度是比较慢的。

不过对于Lambert的BRDF,问题变得简化:它假设了BRDF的结果只取决于入射光,而与出射光无关(即与观察的视角无关),此时相当于表面是一个非常粗糙的完全没有高光的材质。公式:
L o ( p , ω o ) = ∫ Ω c π L i ( p , ω i ) ( n ⋅ ω i ) d ω i = c π ∫ Ω L i ( p , ω i ) ( n ⋅ ω i ) d ω i \begin{aligned} L_o(p,\omega_o) & = \int_\Omega\frac{c}{\pi}L_i(p,\omega_i)(\bf{n} \cdot \omega_i) \it{d} \omega_i \\ & = \frac{c}{\pi}\int_\Omega L_i(p,\omega_i)(\bf{n} \cdot \omega_i) \it{d} \omega_i \\ \end{aligned} Lo​(p,ωo​)​=∫Ω​πc​Li​(p,ωi​)(n⋅ωi​)dωi​=πc​∫Ω​Li​(p,ωi​)(n⋅ωi​)dωi​​

(公式里为什么要除以 π \pi π还不理解,待后续学习)

这样,就可以在所有 n \bf{n} n方向上进行积分,将结果存到一张CubeMap上了,这也是本篇要做的事情。

积分

由于原文作者在shader中将在球面坐标上积分,所以下面以球面坐标的形式改写公式:
L o ( p , θ o , ϕ o ) = c π ∫ Φ ∫ Θ L i ( p , θ i , ϕ i ) c o s ( θ i ) s i n ( θ i ) d θ i d ϕ i L_o(p,\theta_o,\phi_o)=\frac{c}{\pi} \int_\Phi \int_\Theta L_i(p,\theta_i,\phi_i)cos(\theta_i)sin(\theta_i)d\theta_id\phi_i Lo​(p,θo​,ϕo​)=πc​∫Φ​∫Θ​Li​(p,θi​,ϕi​)cos(θi​)sin(θi​)dθi​dϕi​

求解它的一种方式是使用蒙特卡洛,但原文作者选择使用离散化的方式:
L o ( p , θ o , ϕ o ) = c π 2 π N 1 π 2 N 2 ∑ N 1 ∑ N 2 L i ( p , θ i , ϕ i ) c o s ( θ i ) s i n ( θ i ) = π c N 1 N 2 ∑ N 1 ∑ N 2 L i ( p , θ i , ϕ i ) c o s ( θ i ) s i n ( θ i ) \begin{aligned} L_o(p,\theta_o,\phi_o) & = \frac{c}{\pi} \frac{2\pi}{N_1} \frac{\pi}{2N_2} \sum^{N_1} \sum^{N_2} L_i(p,\theta_i,\phi_i)cos(\theta_i)sin(\theta_i)\\ & = \frac{\pi c}{N_1N_2} \sum^{N_1} \sum^{N_2} L_i(p,\theta_i,\phi_i)cos(\theta_i)sin(\theta_i)\\ \end{aligned} Lo​(p,θo​,ϕo​)​=πc​N1​2π​2N2​π​∑N1​​∑N2​​Li​(p,θi​,ϕi​)cos(θi​)sin(θi​)=N1​N2​πc​∑N1​​∑N2​​Li​(p,θi​,ϕi​)cos(θi​)sin(θi​)​

这个公式将原先的积分变成了离散的求和。其中 N 1 N_1 N1​是方位角划分的步数, N 2 N_2 N2​是极角划分的步数

代码

对应的着色器代码如下(原文是hlsl,而下面是glsl):

#version 450
#extension GL_ARB_separate_shader_objects : enable

layout(location = 0) in vec2 InterpolatedPosition;

//UniformBuffer
layout(binding = 0) uniform UniformBufferObject {
	int cubeFace;
} ubo;

//CubeMap纹理
layout(binding = 9) uniform samplerCube ourCubeMap;

out vec4 fColor;

void main()
{
	float PI = 3.1415926;

	//这个函数用来计算“特定表面法线方向”上的出射辐射率,结果将存在CubeMap的对应位置上
	
	//cubeFace代表的是CubeMap中哪一个面
	int cubeFace = ubo.cubeFace;
	//InterpolatedPosition代表“特定表面法线方向”对应在这一个CubeMap的面上的位置,范围(-1~1)

	//下面根据cubeFace来换算normal:“特定表面法线方向”在世界空间下的方向
    vec3 normal = normalize( vec3(InterpolatedPosition.xy, 1) );
    if(cubeFace==2)
        normal = normalize( vec3(InterpolatedPosition.x,  1, -InterpolatedPosition.y) );
    else if(cubeFace==3)
        normal = normalize( vec3(InterpolatedPosition.x, -1,  InterpolatedPosition.y) );
    else if(cubeFace==0)
        normal = normalize( vec3(  1, InterpolatedPosition.y,-InterpolatedPosition.x) );
    else if(cubeFace==1)
        normal = normalize( vec3( -1, InterpolatedPosition.y, InterpolatedPosition.x) );
    else if(cubeFace==5)
        normal = normalize( vec3(-InterpolatedPosition.x, InterpolatedPosition.y, -1) );

	//下面计算球面坐标的轴在世界空间下的方向,
	//其中up是天顶方向,right是方位角基准方向
    vec3 up = vec3(0,1,0);
    vec3 right = normalize(cross(up,normal));
    up = cross(normal,right);

    vec3 sampledColour = vec3(0,0,0);	//逐渐累加的采样得到的颜色值,最后会根据采样次数平均化
    float index = 0;					//采样次数,用于最后的平均化
    for(float phi = 0; phi < 6.283; phi += 0.025)		//采样所有的方位角方向
    {
        for(float theta = 0; theta < 1.57; theta += 0.1)//采样顶半球所有的极角方向
        {
            //下面计算sampleVector:本次采样的方向在世界空间下的方向
			vec3 temp = cos(phi) * right + sin(phi) * up;
            vec3 sampleVector = cos(theta) * normal + sin(theta) * temp;
            
			//从环境贴图中采样,并套入本篇上面讨论的公式
			sampledColour += texture( ourCubeMap, sampleVector ).rgb * 
                                      cos(theta) * sin(theta);
            index = index+1;
        }
    }
	
	//最终结果平均化:
    fColor =  vec4(( PI * sampledColour / index), 1 );
}

其中,cubeFace是UniformBuffer,在C++程序中控制从0递增到5。且每次渲染后都保存为图片。(关于UniformBuffer已经在《图形API学习工程(6):创建并使用UniformBuffer》中做了实践,保存图片的功能在《将OpenGL渲染的结果保存为图片》中做了实践)

//六个面保存的图片:
std::vector<std::string> fileNames = { 
    "right.bmp",
    "left.bmp",
    "top.bmp",
    "bottom.bmp",
    "front.bmp",
    "back.bmp"
     };

//渲染六个面:
for (int i = 0; i < 6; i++)
{
    //更新cubeFace的UniformBuffer:
    UBData.cubeFace = i;
    UpdateUniformBuffer();

    //显示(调用 draw call)
    display();      

    //将当前渲染结果保存为图片
    ScreenShot(fileNames[i]);

    //交换前后缓冲
    glfwSwapBuffers(window);

    //轮询并处理事件
    glfwPollEvents();
}

完整工程代码见 https://codechina.csdn.net/u013412391/cubemaplambertdiffuse

结果

运行程序后,将会依次输出6张贴图。与原始的CubeMap对比如下:
从CubeMap生成Lambert下的环境辐射率图


我用原文中的图片做测试,并仿照他的样子拼合,结果如下:
从CubeMap生成Lambert下的环境辐射率图
但是原文的结果如下:
从CubeMap生成Lambert下的环境辐射率图
目前我并不清楚为什么和原文有差异,比较怀疑的一个原因是原文的CubeMap是HDR,但我这里只是保存了一个低分辨率的图片。

上一篇:凸优化(一)绪论与凸集


下一篇:02梯度下降算法