微表面材质Microfacet BRDF+代码实现+镜面反射+常见问题解决
注意：本文基于games101 作业7中的代码框架
- 如何实现微表面材质？
- 如何提高渲染速度？
- 如何解决白色黑色噪点太多的问题？
- 如何实现镜面反射？
- 如何理解ior在镜面反射的影响
图形学中什么是材质？到底什么是材质？其实就是光线在表面上的各种样式的反射，那什么定义反射呢？当然就是BRDF了，归根到底 材质==BRDF
材质总结
Diffuse
我们现在换一个角度来理解漫反射的表面是怎么反射光的
我们知道，对于一个材质，我们需要知道从任意观测角看这个材质的一个点，这个方向上的颜色是什么，由于漫反射光一般是完全均匀的散射出去，所以从观测方向的w0看到的能量大小：
- fr是brdf，由于是均匀的，所以就是1/π,如果有颜色，那 fr 就是RGB三通道的值Kd比上π
- Li(wi)就是从wi方向来的光能量
- 后面的cos就是处理斜着照的情况
- 积分积的是上半球的每一个入射光的能量
Glossy
这种glossy材质就是反射，但是反射的不彻底，有一点漫反射的味道，这种材质的brdf比较难表示出来，因为不是均匀的反射，所以brdf不是常数1/π
反射折射
需要考虑的有折射率，菲涅尔项 ，全反射。由于反射我们有全向分布函数了，那么折射其实也是有的，叫BTDF，和BRDF 合称BSDF
菲涅尔项
这时一个新的概念，什么是菲涅尔项，其实说明了这么一件事：
-
就是对于一个绝缘体平面，反射的量的多少其实是和你观测方向是有关系的，你垂直观测，其实大部分光都会直接投射折射（前提这个物体透光）了，如果你的观测方向和法线角度很大，那大部分光其实是直接反射的，所以菲涅尔项就告诉你这个比例的值是多少
-
对于导体金属，大部分都是反射
所以我们代码里面实现提供菲涅尔项的函数其实是对于绝缘体的，需要返回一个数值，告诉多少比例是反射的，这怎么计算呢？既然是和表面作用的结果，那计算一定是需要使用到法线N ， 两种介质的折射率n1 n2，光线角度θ，最后我们拟合出一个方程，计算菲涅尔项的值：
其中R0是基准反射率， 就是即使你垂直看下去，还是有一点点反射的，不是全投过去，它的值和介质折射率有关为：
其实真正更准确的计算方法是下面的公式 结果是Reff：
Microfacet BRDF
原理
微表面模型理论中，我们始终假设大部分的材质都是很多小的方向不一的面构成，这些面全是理想的镜面反射的面，也都有自己的法线，朝向不一，就是因为这些小面的朝向的多样性，形成了光线的各种各样的反射，从而诞生了多样的材质
我们从微表面模型来理解之前遇到的材质
Glossy材质
就是有反射，但是反射的范围是扩散的，其实就是说明，微小的表面的法线大致是向上的
Diffuse材质
就是微表面的法线方向不一，很乱，宏观上体现材质是粗糙的
那微表面材质的BRDF 是怎么算的呢？ ， 就是给你一个入射方向，对于其他出射方向的光的能量比例是多少呢？与什么量有关呢？公式如下：（i 观测方向 ， o 入射方向 ， h是半程向量）
解释：
- F（i , h ）菲涅尔项：输入的是观测方向和半程向量，返回一个比例，告诉反射折射各自的占比
- D（h）法线分布函数：反应了微表面的法线方向的分布函数，告诉你所有微表面有多少比例是可以反射你这个光的，所以需要传入半程向量，在D里面计算有多少比例的微表面的法线是大致和你这个半程向量重合的，最后返回这个比例
- G(i , o , h)几何衰减项:处理微表面之间自己遮挡的情况，就是我在D里算出来有这么多的能量会反射，但是由于微表面之间自己遮挡，会少一点，G就是处理这件事的
微表面模型有很多种，但是都基于小平面的这一套
BRDF的计算
那上面的F菲涅尔项 ， G几何衰减项 ， D 法线分布函数 的计算公式是什么呢？
F菲涅尔项
菲涅尔项（观测方向 ， 法向 ， 折射率 ），需要计算出：在没有发生全反射的情况下，反射和折射能量各自的占比是多少，，这里上文已经给出计算公式
其中R0是基准反射率， 就是即使你垂直看下去，还是有一点点反射的，不是全投过去，它的值和介质折射率有关为（空气的 ior 一般我们考虑成1 ）：
这里直接使用作业7里面的，I 是光来的方向（注意作业7中间的入射方向和课程讲的是 相反的，所以，注意正负号的问题！），N是宏观的法线 ， ior是折射率 ， Kr就是需要我们需要计算的反射和折射的比率
void fresnel(const Vector3f &I, const Vector3f &N, const float &ior, float &kr) const
{
float cosi = clamp(-1, 1, dotProduct(I, N));
float etai = 1, etat = ior;
if (cosi > 0) { std::swap(etai, etat); }
// Compute sini using Snell's law
float sint = etai / etat * sqrtf(std::max(0.f, 1 - cosi * cosi));
// Total internal reflection
if (sint >= 1) {
kr = 1;
}
else {
float cost = sqrtf(std::max(0.f, 1 - sint * sint));
cosi = fabsf(cosi);
float Rs = ((etat * cosi) - (etai * cost)) / ((etat * cosi) + (etai * cost));
float Rp = ((etai * cosi) - (etat * cost)) / ((etai * cosi) + (etat * cost));
kr = (Rs * Rs + Rp * Rp) / 2;
}
// As a consequence of the conservation of energy, transmittance is given by:
// kt = 1 - kr;
}
G几何衰减项
有很多中集合衰减的函数G，其中分离遮蔽阴影最常用，这种方法和后面提到的法线分布函数GGX共同组成是目前UE4中正在使用的方案
其中 i 是观测方向 ， o是入射方向，且Gs如下：
其中
注意α表示的是粗糙度，小说明反光，大说明粗糙！！
代码实现集合衰减项，i 是入射，o是出射 ， a 是粗糙度，一般设置0 - 20 之间，注意，只要有点乘，就需要约束在0 到1 之间！！！
float G_s(float NdotV, float a)
{
float k = (a + 1) * (a + 1) / 8;
return NdotV / (NdotV* (1- k) + k );
}
float G(Vector3f i, Vector3f o, Vector3f N , float a)
{
float NdotI = clamp(0.0, 1.0, dotProduct(N, i));
float NdotO = clamp(0.0, 1.0, dotProduct(N, o));
return G_s(NdotI, a) * G_s(NdotO , a);
}
D 法线分布函数
业界较为主流的法线分布函数是GGX
α表示的是粗糙度
代码：r 是粗糙度
float D_GGX(Vector3f h, float r, Vector3f N)
{
float r2 = r * r;
float NdotH = clamp(0.0, 1.0, dotProduct(N, H)); //有点乘就需要约束！！
float NdotH2 = NdotH * NdotH;
float res = r2 / (M_PI * (NdotH2 * (r2 - 1) + 1) * (NdotH2 * (r2 - 1) + 1));
return res;
}
作业7中部分函数变化
这样我们的材质中计算BRDF的函数，从原来的只有漫反射到现在实现了微表面材质的了，作业7中Material 里面的成员函数eval 需要改
同样需要注意：
由于作业7的入射光我吗是从外设向着色点的，所以涉及到入射光wi的就需要加上一个负号！不然全黑，如果物体全白就是G里面的wi忘记写负号了！
点乘需要约束范围 ，记得加一个新的枚举类型 ， 防止分母除以0 需要加一个0.0001 ， 同时取绝对值
//wi是光来的方向 //wo是光去的方向
Vector3f Material::eval(const Vector3f &wi, const Vector3f &wo, const Vector3f &N)//其实这就是bsdf
{
switch(m_type){
case DIFFUSE:
{
float cosalpha = dotProduct(N, wo); //wo是光去的方向（我们光源是从像素出发的
if (cosalpha > 0.0f) {
Vector3f diffuse = Kd / M_PI; //之前全是diffuse材质，直接返回的是Kd/π
return diffuse;
}
else
return Vector3f(0.0f);
break;
}
case MIRCO:
{
float cosalpha = dotProduct(N, wo); //wo是观测方向
if (cosalpha > 0.0f) {
Vector3f h = (-wi + wo).normalized();
float F = 0.f;
fresnel(-wi, N, ior, F);
float down = 4 * fabs(clamp(0.0, 1.0, dotProduct(N, -wi)) * clamp(0.0, 1.0, dotProduct(N, -wi))) + 0.00001;
float up = F * G( -wi , wo , N , roughness) * D_GGX( h , roughness , N); //F * G * D
return (up / down) * Kd; //加上Kd就实现了颜色，当然主函数记得给材质设置Kd
}
else
return Vector3f(0.0f);
break;
}
}
}
这样返回微表面材质下的BRDF的函数也完成了，这样在光追渲染时，如果这个材质是微表面，就会计算对应的brdf了
注意主函数需要自己建立几个球，写好材质，我自己加了两个球，上了色（Kd）
Material* Microfacet1 = new Material(MIRCO, Vector3f(0.0f));
Microfacet1->ior = 5;
Microfacet1->roughness = 0.06;
Microfacet1 ->Kd = Vector3f(0.14f, 0.60f, 0.091f);
Material* Microfacet2 = new Material(MIRCO, Vector3f(0.0f));
Microfacet2->ior = 5;
Microfacet2->roughness = 0.2;
Microfacet2 ->Kd = Vector3f(0.67f, 0.065f, 0.05f);
Sphere b1(Vector3f(160, 120, 150), 120, Microfacet1);
Sphere b2(Vector3f(400 , 120, 350), 120, Microfacet2);
scene.Add(&floor);
//scene.Add(&shortbox);
scene.Add(&b1);
scene.Add(&b2);
//scene.Add(&tallbox);
scene.Add(&left);
scene.Add(&right);
scene.Add(&light_);
注意球的类中的求交点的函数有浮点数精度问题，需要改一点 getIntersection 函数！！，
Intersection getIntersection(Ray ray)
{
Intersection result;
result.happened = false;
Vector3f L = ray.origin - center;
float a = dotProduct(ray.direction, ray.direction);
float b = 2 * dotProduct(ray.direction, L);
float c = dotProduct(L, L) - radius2;
float t0, t1;
if (!solveQuadratic(a, b, c, t0, t1)) return result;
if (t0 < 0) t0 = t1;
if (t0 < 0) return result;
if (t0 > 0.5) //这里这里！！！！
{
result.happened = true;
result.coords = Vector3f(ray.origin + ray.direction * t0);
result.normal = normalize(Vector3f(result.coords - center));
result.m = this->m;
result.obj = this;
result.distance = t0;
}
return result;
}
效果图
像下面这样的，如果你的渲染图有很多白色黑色的噪点之类的，还有边缘发光发亮的，估计是在castRay 里面没有限制颜色的范围倒置的！
需要在castRay 的结尾使用一下手段约束一下颜色的范围
Vector3f res = L_dir + L_indir;
res.x = clamp(0.0, 1.0, res.x);
res.y = clamp(0.0, 1.0, res.y);
res.z = clamp(0.0, 1.0, res.z);
return res;
这样就不会有白色黑色的噪点了，图片质量大福提升！！！图的spp只有18呢
**大幅提高渲染速度 **
除了使用多线程，在随机数的生成上面 ，需要将函数和随机数的构造过程全部设置为static的，这样只会构造一次，而不是每次需要随机数的时候就构造一次，不然相当的耗时！！
inline static float get_random_float()
{
static std::random_device dev;
static std::mt19937 rng(dev());
static std::uniform_real_distribution<float> dist(0.f, 1.f); // distribution in range [1, 6]
return dist(rng);
}
实现反射
现在需要实现这种效果，没错就是镜面反射
这时候，由于需要实现镜面反射，所以，在光打到表面之后，反射的方向的采样就必须全部是反射方向，所以，Material::sample需要这样改
Vector3f Material::sample(const Vector3f &wi, const Vector3f &N){ //上半球的随机采样
switch(m_type){
case MIRCO:
case DIFFUSE:
{
// uniform sample on the hemisphere
float x_1 = get_random_float(), x_2 = get_random_float();
float z = std::fabs(1.0f - 2.0f * x_1);
float r = std::sqrt(1.0f - z * z), phi = 2 * M_PI * x_2;
Vector3f localRay(r*std::cos(phi), r*std::sin(phi), z);
return toWorld(localRay, N);
break;
}
case REFLC: //这里这里！！
{
Vector3f localRay = reflect(wi, N).normalized();
return localRay;
break;
}
}
}
同时PDF也需要改的，由于是镜面反射，所以反射至采样方向的概率是1
float Material::pdf(const Vector3f &wi, const Vector3f &wo, const Vector3f &N){
switch(m_type){
case MIRCO:
case DIFFUSE:
{
// uniform sample probability 1 / (2 * PI)
if (dotProduct(wo, N) > 0.0f)
return 0.5f / M_PI;
else
return 0.0f;
break;
}
case REFLC:
{
if (dotProduct(wo, N) > 0.0001f)
return 1.0;
return 0.0;
beak;
}
}
}
那么镜面反射的BRDF怎么计算呢？
由于影响镜面反射的强度只是和角度和入射光线共同组成的菲涅尔项有关 ， 为了保证最终结果只和菲涅尔项和反射光线有关，brdf里要抵消掉cosθ的影响，所以
brdf = 1 / cosθ
乘上由于菲涅尔项造成的比率 K
*brdf = K
如果有颜色的镜子， 就是乘上Kd
*brdf = Kd
于是eval如下
case REFLC:
{
float cosa = dotProduct(N, wo);
if (cosa > 0.0001f)
{
float K = 0.0;
fresnel(wi, N ,ior , K);
Vector3f t = (1.0 / cosa) ;
t *= Kd;
return t * K;
}
break;
}
同时， 需要在castRay 改一下，渲染着色需要和镜面反射的区分开来，镜面反射就不需要直接对光源采样了，不然可以被光照的地方镜子是全白的！
Vector3f Scene::castRay(const Ray &ray, int depth) const
{
// TO DO Implement Path Tracing Algorithm here
Intersection hit = intersect(ray);
Vector3f L_dir, L_indir;
if (!hit.happened)
return L_dir;
if (hit.m->hasEmission())return hit.m->getEmission();
//接下来说明打到的是物体
//将交点的信息取出来
Vector3f p = hit.coords;
Material* m = hit.m;
Vector3f N = hit.normal;
Vector3f w0 = ray.direction; //注意这个入射方向是光射向着色点的！！！！！！
w0.normalized();
switch (m->getType())
{
case REFLC:
{
float prr = get_random_float();
if (prr > RussianRoulette) return L_indir;
Vector3f wi = m->sample(w0, N).normalized();
Ray r(p + 0.0001, wi);
Intersection hit2 = intersect(r);
if (hit2.happened)
{
Vector3f brdf2 = m->eval(w0, wi, N);
float cos_theta3 = dotProduct(wi, N);
float pdf_ = m->pdf(w0, wi, N);
if (pdf_ > 0.0001)
{
L_indir = castRay(r, depth + 1) * brdf2 * cos_theta3 / pdf_ / RussianRoulette;
}
}
break;
}
case DIFFUSE:
case MIRCO:
{
float pdf_L;
Intersection light_hit;
sampleLight(light_hit, pdf_L);//随机在所有光源上选出一个点，并且得到这个点的pdf
//看看这个点可不可以照射到摄像头发出的光线打到的物体的点hit上
Vector3f sa_l_coord = light_hit.coords;
Vector3f sa_l_N = light_hit.normal;
Vector3f dir_hit_to_sa_l = (sa_l_coord - p).normalized();
Vector3f intencity_of_L = light_hit.emit;
float d = (sa_l_coord - p).norm();
Ray isBlock(p + 0.001, dir_hit_to_sa_l);
Intersection hi = intersect(isBlock);
if (hi.happened && hi.m->hasEmission()) {
Vector3f brdf = m->eval(w0, dir_hit_to_sa_l, N);
float cos_theta1 = dotProduct(N, dir_hit_to_sa_l);
float cos_theta2 = dotProduct(sa_l_N, -dir_hit_to_sa_l);
L_dir = (intencity_of_L * brdf * cos_theta1 * cos_theta2 / std::pow(d, 2)) / pdf_L;
}
//计算间接光照
float prr = get_random_float();
if (prr < RussianRoulette)
{
Vector3f wi = m->sample(w0, N).normalized();
Ray ri(p, wi);
Intersection hit2 = intersect(ri);
if (hit2.happened && hit2.m->hasEmission() == false)
{
Vector3f brdf2 = m->eval(w0, wi, N);
float cos_theta3 = dotProduct(wi, N);
float pdf_ = m->pdf(w0, wi, N);
L_indir = castRay(ri, depth + 1) * brdf2 * cos_theta3 / pdf_ / RussianRoulette;
}
}
}
}
Vector3f res = L_dir + L_indir;
res.x = clamp(0.0, 1.0, res.x);
res.y = clamp(0.0, 1.0, res.y);
res.z = clamp(0.0, 1.0, res.z);
return res;
}
接下来可以渲染了
效果如下，有Kd颜色镜面反射的：
没颜色Kd的镜面反射
如何理解 ior 对镜面反射的影响
ior是折射率 ，其实与其说我们对球体实现的是镜面反射，不如说是实现了一个类似水的球体，有折射有反射，只是由于折射我们没有实现，所以球里面出来的光是黑的，我们看到的就只有球的反射的样子
ios越大，由于**菲涅尔项 ， **反射比率越高， 镜子越亮 ； ios越小，由于大多光线透射了，而且我们作业7没有实现折射，所以从里面出不来光，所以 ios越小 ， 球越黑
就这样我们实现了微表面材质，实现了镜面反射，喜欢的收藏点赞关注吧
