( Precomputed Radiance Transfer ) PRT

PRT次表面散射

次表面散射(SubSurface Scattering)是指光线照射到物体后，进入物体内部，经过在物体内的散射从物体表面的其他顶点/像素离开物体的现象。在透明/半透明的材质上，这种效果比较明显，比如说是皮肤，玉等等。和BRDF的定义一样我们将存在次表面散射情况下的反射光亮度/照度定义为 BSSRDF。我们知道BRDF是BSSRDF的一种特殊情况，即所有的光线从一个顶点/像素进入物体内部也由该顶点/像素离开物体。从上面次表面散射的定义我们知道我们考虑该种情况下的反射不仅要考虑该顶点/像素对于入射光线的反射，还需要考虑表面其他的顶点/像素的对于该顶点/像素的散射光线。而BSSRDF是一个全局的考虑，它包含了BRDF的情况。

包含次表面散射的情况下，渲染公式是

Lo =∫∫BSSRDF*Li(xi,wi)*cosθi*dwi*dA(xi)

这里BSSRDF是xi,xo,wi,wo的函数。我们简化这个函数，使得它和wo无关，这样在投影阶段我们就可以使用Vector的系数来表示。也就是我们对于片面A上的每个顶点/像素计算对于顶点/像素xo的散射贡献度。我们首先使用蒙特卡罗积分法对于最外面的对于片面A的积分，我们在片面A上随即采样N，计算每次的散射值。p(xi)=1/S(A)。得到

Lo ≈(1/N)*∑(∫BSSRDF*Li(xi,wi)*cosθi*dwi)/p(xi)
≈(S(A)/N)*∑∫BSSRDF*Li(xi,wi)*cosθi*dwi

接下来我们将Li(xi,wi)和BSSRDF* cosθi预先投影到Spherical Harmonic,得到系数Tji和Lji。运行时我们只需要

Lo ≈(A/N)*∑∑Tji*Lji

这里我们需要注意的是对于我们的Precomputed Radiance Transfer计算我们需要简化上面的计算。我们可以假设在片面A上，各处的入射光一致，那么我们就可以使用该顶点/像素的Li来代替片面上其他顶点/像素的Lji.这样我们最终的计算就可以表示为:

Lo ≈(A/N)*∑(Li*∑Tji)

//
//PS:对于DX9中PRT处理SSS我始终找不到相关实现的方法资料。
//相关BSSRDF的计算可以参考
//“A Practical Model for Subsurface Light Transport”中的相关推导。
//
作者Blog：http://blog.csdn.net/sevecol/

PRT计算光照
直接光照和间接光照

我们现在的任务是要在上面直接光照和阴影的基础上加入的是间接光照，也就是需要加上其他的表面反射光线的照射对于我们当前计算顶点/像素的光照效果。我们可以和Path Raytrace所采用的方法一样，首先我们的渲染公式分为两个部分，直接照射和间接照射。

L =L1+L2
=L1+∫BRDF*L(x1)*cosθ*V*dω

这里L1表示直接光照，L2表示间接光照。L(x1)表示其他表面反射的光亮度。

和上面的介绍的两种情况一样，我们假设我们的顶点/像素所在的平面是理想的 Diffuse表面。我们接下来要解决得是我们公式中的间接光照部分。我们对于间接光照部分采用蒙特卡洛积分法来计算，于是可以得到

L =L1+∫BRDF*L(x1)*cosθ*dω
≈L1+(ρ1/PI)*(1/N)*∑(L(x1)*cosθ)/p(x1)

通过上面的表达式，我们可以得出我们计算间接光照部分的步骤是:首先我们以当前需要计算反射的顶点/像素为中心，在顶点/像素和所在表面的法线确定的半球内随机发射 N条光线和场景以及其他物体做碰撞检测。如果没有和场景或者其他物体都没有碰撞，那么这个方向不存在间接光照，也可以表示为L(x1)=0。反之就表示这个方向存在间接光照。然后我们就需要计算出这个交点朝该方向的反射光亮度。而对于这个交点来说朝这个方向的反射光亮度等于这个交点对于直接光照和间接光照的反射光亮度。我们在这里仅仅使用这个交点直接光照的反射光亮度来近似的代替交点朝这个方向的反射光亮度。这样的话我们就需要对于交点计算直接光照反射光亮度。我们使用上面提到的直接光照+阴影的方法来计算交点直接光照的反射光亮度。于是我们进一步得到

L ≈L1+(ρ1/PI)*(1/N)*∑(L(x1)*cosθ)/p(x1)
≈L1+(ρ1/PI)*(1/N)* ∑((∫BRDF*L(x1)*cosθ*dω)* cosθ)/p(x1)
≈L1+(ρ1/PI)*(1/N)* ∑∑((Li*T2i)* cosθ)/p(x1)
≈(ρ/PI)*∑Li*T1i+(ρ1/PI)*(1/N)* ∑∑((Li*T2i)* cosθ)/p(x1)

我们在当前计算顶点/像素和所在平面的法线所确定的半球随机采样，在这样的情况下p(x1)=1/(2*PI),带入上面最后的表达式得到

L ≈(ρ/PI)*∑Li*T1i+(ρ1/PI)*(1/N)*2*PI*∑∑((Li*T2i)* cosθ)
≈(ρ/PI)*∑Li*T1i+(ρ1/PI)*(1/N)*2*PI*∑(Li*∑(T2i* cosθ))

通过上面的推导，我们可以看出我们在计算直接光照和间接光照效果的时候，我们计算顶点/像素的T(x)的系数Ti，我们首先需要计算出顶点/像素在直接光照+阴影效果下的T(x)系数T1i。然后在顶点/像素和所在的平面的法线确定的半球内我们随机采样N条光线和场景以及物体作碰撞检测。对于每个交点，我们使用插值计算出该交点在直接光照效果下的T(x)系数。将我们得到所有的交点的系数乘上该光线和该交点所在的平面的法线的余铉值cosθ，然后求出平均值T2i。我们最终的直接光照和间接光照的T(x)系数了就等于我们上面计算出来的两组系数的和。Ti=T1i+T2i。

最终的渲染还是使用预计算的Spherical Harmonic Lighting和Radiance Transfer的系数之间的做点乘操作得到了最终需要计算的顶点/像素的反射光的亮度值。

作者Blog：http://blog.csdn.net/sevecol/


PRT计算Volume Rendering

我们知道离开Volume时候的光照亮度
Lo=Li*T+∫T*(Le(p,w)+σs(p,w)*∫p(w1,w)*L(w1)*dw1)*dt
其中T=exp(-∫(σs(p,w)+σa(p,w))*dt).
我们注意到公式中的In-Scattering部分
∫σs(p,w)*∫p(w1,w)*L(w1)*dw1)*dt
这里有两次积分过程。第一次是计算四周的光线散射到视线方向，第二次积分是对光线传播途中进行积分。我们知道第一次积分的时候的四周的光线其实也是多次散射的结果，我们这里适当的简化，我们假设我们在
第一次积分的时候只考虑直接照射的光线经过吸收和Out-Scattering后剩余的亮度，而不考虑散射来的光线亮度。这样就可以得到简化后的In-Scattering部分
∫σs(p,w)*∫p(w1,w)*L1(w1)*dw1)*dt //L1(w1)代表w1方向的只考虑吸收和Out-Scattering过程光线亮度
我们知道只考虑吸收和Out-Scattering的时侯 L1(w1)=Ls*Ts(w1),这里Ls表示光源的直接照射，Ts(w1)代表光线到达散射微粒的吸收和Out-Scattering过程。于是In-Scattering部分为
∫σs(p,w)*∫p(w1,w)*Ts(w1)*Ls(w1)*dw1)*dt
我们观察内部的对于球形的积分，我们的积分空间是以微粒为中心所在的球。那么我们对于这个积分可以使用SH来解开成为dot操作。我们将积分号内部分成两个部分,p(w1,w)*Ts(w1)和Ls(w1)。我们假设整个微粒所在的Volume对于环境来说很小，那么每一个微粒的Ls(w1)都是一样的，不同的就是在p(w1,w)*Ts(w1)部分。我们将这两个部分投影到SH上面可以分别得到系数。对于每一个微粒来说他们的Ls(w1)的系数是一样的，p(w1,w)*Ts(w1)的系数是不一样的，而且对于每一个w都有一组不一样的系数。我们可以把对于不同的w的系数保存在cube map中以备查找。这样我们碰到一个问题，就是每个微粒都有cube map,这样的结构对于硬件是不友好的。我们只能继续简化我们的操作。我们将上面提到的Ts(w1)忽略掉，那么对于所有的微粒就可以用同样的cube map来保存了(可能会需要多张，因为精度的需要)。
我们在下来看看最外面的积分，我们首先假设σs(p,w)在传播途中是恒量。我们可以使用蒙特卡罗积分法来近似的计算这个积分，我们取视线穿过Volume的两头的作为我们的采样,那么我们只需要计算出这两点微粒的结果和/2就可以得到近似值了。
通过上面的操作，我们有了对于微粒的N组系数(p(w1,w)投影得到)Ti(w)保存在cube map中，靠w也就是视线的方向来获取。同时还有周围环境的一组Ls(w1)的系数Li。我们通过一些技巧还能获得视线方向的Volume深度d.那么计算的时候我们首先渲染一边Volume到Texture,这次渲染我们把Volume的前面和后面的Ti(w)相加。接下来我们把这个相加得到的值和环境的Li作dot操作。并乘上σs*d/2得到了In-Scattering的结果。
最后把视线方向的吸收,Out-Scattering和上面计算得到的In-Scattering加起来就能得到我们所需的近似值了。
还没有实践不知道最终的效果会是怎么样的。