鉴于国内网络能够搜到关于SVO的文章较少，阅读国外Paper十分费心费力，我准备将最近阅读的 VCTRenderer 代码中所学到的东西记录下来，作为自己的知识储备。

二、VCTRenderer的大体渲染流程

三、阴影贴图的计算

VCTRenderer中为了达到软阴影的效果，没有使用PCF，而是使用了EVSM——Exponential Variance Shadow Map，该算法结合了Exponential Shadow Map和Variance Shadow Map，能够在达到软阴影效果的同时有效减少光渗透现象。

1、Exponential Shadow Map

该公式中， d(x) 代表的是在Fragment Shader中将 x 变换到光源空间时所得到的深度值， z(p) 代表的是从阴影贴图中纹理坐标为 p 的位置采样得到的深度值。 f 则是映射关系，比如硬阴影的映射关系如下：

if (d(x) < z(p))
    return 1;//被遮挡
    return 0;//未被遮挡

为了得到软阴影，我所知道的一般有两种处理方式，第一种是不对阴影贴图做任何处理，在Fragment Shader中通过PCF多次采样进行模糊处理，代码如下：

shadow = 0.0;
vec2 texelSize = 1.0 / textureSize(ShadowMap, 0);
for(int x = -1; x <= 1; ++x)
    for(int y = -1; y <= 1; ++y)
        float pcfDepth = texture(ShadowMap, projCoords.xy + vec2(x, y) * texelSize).r; 
        shadow += currentDepth - bias > pcfDepth  ? 1.0 : 0.0;        
shadow /= 9.0;

第二种方式是直接对生成后的阴影贴图进行模糊处理，例如高斯模糊，然后在Fragment Shader中进行采样，但是由于进行了模糊处理——或者说卷积处理：

那么将 f(d(x),z(p)) 替换 g(p-q) 得到下式：

因为我们只需要对 z(p) 进行卷及处理，为了将 d(x) 与 z(p) 分离，Exponential Shadow Map使用了指数函数代替了Shadow Map的数学模型：

将上式替换 f(d(x),z) ，且由于我们只关心Shadow Map可以认为 d(x) 是一个常数（因为它在Fragment Shader中对于同一个像素来说，值恒定），可以得到下式：

由此可见，成功将 d(x) 与 z(p) 分离，我们接下来只需要对 z(p) 进行卷积操作即可。

PS：该节内容主要参考自 该博客 。

2、Variance Shadow Map

3、Exponential Variance Shadow Map

Exponential Variance Shadow Map是Exponential Shadow Map与Variance Shadow Map相结合的产物，网络上的相关资料比较少（没有找到），这里直接贴上VCTRenderer中的实现方式：

depth_texture.frag中：

vec2 WarpDepth(float depth)
    // rescale depth into [-1, 1]
    depth = 2.0 * depth - 1.0;
    float pos =  exp( exponents.x * depth);//exponents的值为vec2(40,5)
    float neg = -exp(-exponents.y * depth);
    return vec2(pos, neg);
vec4 ShadowDepthToEVSM(float depth)
	vec2 moment1 = WarpDepth(depth);
	vec2 moment2 = moment1 * moment1;
	return vec4(moment1, moment2);
void main()
    vec4 diffuseColor = texture(diffuseMap, texCoord);
    if (diffuseColor.a <= alphaCutoff) { discard; }
	outColor = ShadowDepthToEVSM(gl_FragCoord.z);
计算光照的Shader中：
 
float linstep(float low, float high, float value)
    return clamp((value - low) / (high - low), 0.0f, 1.0f);
float ReduceLightBleeding(float pMax, float Amount)  
    return linstep(Amount, 1, pMax);  
vec2 WarpDepth(float depth)
    depth = 2.0f * depth - 1.0f;
    float pos = exp(exponents.x * depth);
    float neg = -exp(-exponents.y * depth);
    return vec2(pos, neg);
float Chebyshev(vec2 moments, float mean, float minVariance)
    if(mean <= moments.x)
        return 1.0f;
        float variance = moments.y - (moments.x * moments.x);
        variance = max(variance, minVariance);
        float d = mean - moments.x;
        float lit = variance / (variance + (d * d));
        return ReduceLightBleeding(lit, lightBleedingReduction);
// 指数方差阴影贴图计算可见度
float Visibility(vec3 position)
    //变换到光源视角空间
    vec4 lsPos = lightViewProjection * vec4(position, 1.0f);
    // 避免除以0造成算数错误
    if(lsPos.w == 0.0f) 
	return 1.0f;
    // 转换到NDC空间
    lsPos /= lsPos.w;
    vec4 moments = texture(shadowMap, lsPos.xy);
    // 偏移Z值避免发生light acen现象
    vec2 wDepth = WarpDepth(lsPos.z - 0.0001f);
    // 计算一个可接受的最小方差
    vec2 depthScale = 0.0001f * exponents * wDepth;
    vec2 minVariance = depthScale * depthScale;
    // 分别对正数项和负数项使用切比雪夫不等式计算可见度
    float positive = Chebyshev(moments.xz, wDepth.x, minVariance.x);
    float negative = Chebyshev(moments.yw, wDepth.y, minVariance.y);
    // 选取最小的可见度
    return min(positive, negative);
四、模型体素化
 
1、什么是体素化
 
如上图所示，用一个个小立方体将模型表示出来即为体素化过程，比如《我的世界》这款游戏就可以当做是用体素将一个个模型构造出来。
 
2、使用正交投影矩阵，以及判断投影方向
 
因为体素本质上是一个个轴对齐的立方体，为了将模型的每个三角面用一堆体素表示出来，我们需要使用正交投影，如下图：
 
但我们都清楚的是，一个立方体有六个面，难道我们需要进行六次投影吗？答案是否定的，首先由于是正交投影，因此投影到上面和下面的结果是一样的，左右面、前后面同理，因此可以首先排除三个面，再有就是裂缝的问题，如下两张图：
 
第一幅图是投影到右面得到的结果，第二幅图是投影到上面得到的结果，很明显的看到第二幅图体素化的结果中，体素与体素之间有明显的裂缝，为了解决这个问题，我们要让模型投影后的面积尽可能大，解决方法很简单，只需要判断顶点的法线的最大分量即可，代码如下：
 
int CalculateAxis()
	vec3 p1 = gl_in[1].gl_Position.xyz - gl_in[0].gl_Position.xyz;
	vec3 p2 = gl_in[2].gl_Position.xyz - gl_in[0].gl_Position.xyz;
	vec3 faceNormal = cross(p1, p2);
	float nDX = abs(faceNormal.x);
	float nDY = abs(faceNormal.y);
	float nDZ = abs(faceNormal.z);
	if( nDX > nDY && nDX > nDZ )
	    return 0;//投影到X轴方向上
	else if( nDY > nDX && nDY > nDZ  )
	    return 1;//投影到Y轴方向上
	    return 2;//投影到Z轴方向上
3、解决孔洞问题
 
选择了投影方向后，如果我们直接进行投影操作会出现孔洞问题，其原因是我们最终要在Fragment Shader中完成体素化过程，而在执行Fragment Shader之前会有一个光栅化过程，该过程中会判断模型是否覆盖了像素的中心位置，如果覆盖了则表示该像素需要执行后续的Fragment Shader，否则则不需要，因此模型的一些部分可能会被丢失，为了解决这个问题，需要使用一种称之为保守光栅化的算法。
 
如上图所示，第一步我们需要给三角形面片构造一个包围盒，由于这是投影后的三角形面片，所以只需要一个二维AABB包围盒即可：
 
vec4 AxisAlignedBoundingBox(vec4 pos[3], vec2 pixelDiagonal)
	vec4 aabb;
	aabb.xy = min(pos[2].xy, min(pos[1].xy, pos[0].xy));
	aabb.zw = max(pos[2].xy, max(pos[1].xy, pos[0].xy));
	// 让包围盒稍微大一点
	aabb.xy -= pixelDiagonal;
	aabb.zw += pixelDiagonal;
	return aabb;
第二步，计算三角面片投影后所在的平面方程的系数
 
vec4 trianglePlane;
trianglePlane.xyz = cross(pos[1].xyz - pos[0].xyz, pos[2].xyz - pos[0].xyz);
trianglePlane.xyz = normalize(trianglePlane.xyz);         //计算平面的法线
trianglePlane.w = -dot(pos[0].xyz, trianglePlane.xyz);    //计算平面距原点的距离
 
第三步，计算向外扩展后的三角形面片的顶点坐标
 
由于是正交投影，我们其实可以先忽略顶点坐标的Z分量，只关心X、Y分量，并以三角形面片的每条边的两个端点以及原点构造出三个平面。
 
已知平面方程为：
 
忽略Z分量，并且该平面过原点，可知：
 
化简后可得：
 
因此我们首先构造出三个过原点的平面：
 
vec3 planes[3];
planes[0] = cross(pos[0].xyw - pos[2].xyw, pos[2].xyw);
planes[1] = cross(pos[1].xyw - pos[0].xyw, pos[0].xyw);
planes[2] = cross(pos[2].xyw - pos[1].xyw, pos[1].xyw);
planes[0].z -= dot(halfPixel, abs(planes[0].xy));
planes[1].z -= dot(halfPixel, abs(planes[1].xy));
planes[2].z -= dot(halfPixel, abs(planes[2].xy));
 
此时planes[i].x即为公式中的a，planes[i].y即为公式中的b，planes[i].z即为公式中的d。
 
接下来，求出这三个平面相交所得到的交线，并除以W分量得到Xndc与Yndc：
 
vec3 intersection[3];
intersection[0] = cross(planes[0], planes[1]);
intersection[1] = cross(planes[1], planes[2]);
intersection[2] = cross(planes[2], planes[0]);
intersection[0] /= intersection[0].z;
intersection[1] /= intersection[1].z;
intersection[2] /= intersection[2].z;
 
最后带入到第二步求得的三角形面片所在的平面公式中算出Z分量，即为：
 
其中，intersection[i].x即为公式中的Xndc，intersection[i].y即为公式中的Yndc，trianglePlane[i].x即为公式中的a，trianglePlane[i].y即为公式中的b，trianglePlane[i].z即为公式中的c，trianglePlane[i].w即为公式中的d，Wndc为1。
 
所以有以下代码：
 
float z[3];
z[0] = -(intersection[0].x * trianglePlane.x + intersection[0].y * trianglePlane.y + trianglePlane.w) / trianglePlane.z;
z[1] = -(intersection[1].x * trianglePlane.x + intersection[1].y * trianglePlane.y + trianglePlane.w) / trianglePlane.z;
z[2] = -(intersection[2].x * trianglePlane.x + intersection[2].y * trianglePlane.y + trianglePlane.w) / trianglePlane.z;
pos[0].xyz = vec3(intersection[0].xy, z[0]);
pos[1].xyz = vec3(intersection[1].xy, z[1]);
pos[2].xyz = vec3(intersection[2].xy, z[2]);
 
4、根据保守光栅化后的顶点位置计算对应的体素坐标
 
//从ndc空间变换到世界空间
vec4 voxelPos = vec4(pos[i].xyz * pos[i].w, pos[i].w);
voxelPos.xyz = (viewProjectionI * voxelPos).xyz;
//因为imageStore这个3D纹理的大小是volumeDimension^3（比如256*256*256），
//所以要先将世界坐标转换到模型的AABB包围盒的坐标系当中
//然后除以voxelSize，获得该坐标对应的是哪一个体素，voxelSize=voxelScale * volumeDimension
voxelPos.xyz -= worldMinPoint;
voxelPos *= voxelScale;
Out.wsPosition = voxelPos.xyz * volumeDimension;
 
5、在Fragment Shader中的一些处理
 
经过上述四个步骤后，已经做到将模型分割成一个个小立方体，由于我们使用3D纹理贴图来存储每个体素的信息，所以需要将每个立方体的信息进行处理（颜色、法线、金属度、粗糙度）并通过imageStore函数存储到3D贴图当中。
 
要是用的3D纹理如下（可以添加其他的3D纹理来存储其他的信息）：
 
//volatile保证每次采样该变量时都要从显存中重新读取而不是使用寄存器中的备份
//coherent保证每次采样数据的时候不会从cache中读取数据而是从内存中读取，维护cache与内存的数据一致性
//综上所述volatile coherent保证着色器程序能够在每次采样时获得最新的数据
//题外话memory barrier——内存屏障，对于多线程操作，使用该函数可以保证访问内存的顺序
layout(binding = 0, r32ui) uniform volatile coherent uimage3D voxelAlbedo;
layout(binding = 1, r32ui) uniform volatile coherent uimage3D voxelNormal;
layout(binding = 2, r32ui) uniform volatile coherent uimage3D voxelEmission;
layout(binding = 3, r8) uniform image3D staticVoxelFlag;
 
首先，由3中的图片，我们知道扩展后的三角形面片，有一部分是在AABB包围盒之外的，对此我们要将其裁减掉：
 
if( In.position.x < In.triangleAABB.x || In.position.y < In.triangleAABB.y || In.position.x > In.triangleAABB.z || In.position.y > In.triangleAABB.w )
	discard;
接下来是非常重要的一步，因为同一个体素往往对应多个不同的三角形面片，因此我们需要求得平均值：
 
void imageAtomicRGBA8Avg(layout(r32ui) volatile coherent uimage3D grid, ivec3 coords, vec4 value)
    value.rgb *= 255.0;                 
    uint newVal = convVec4ToRGBA8(value);
    uint prevStoredVal = 0;
    uint curStoredVal;
    uint numIterations = 0;
    //imageAtomicCompSwap：从grid中纹理坐标为coords处采样，若与prevStoredVal相等则将值替换成newVal，返回值是未替换之前的纹理坐标所对应的值
    //这部分代码的作用是当发生重叠的时候要计算出一个平均值
    while((curStoredVal = imageAtomicCompSwap(grid, coords, prevStoredVal, newVal)) != prevStoredVal && numIterations < 255)
        prevStoredVal = curStoredVal;
        vec4 rval = convRGBA8ToVec4(curStoredVal);
        rval.rgb = (rval.rgb * rval.a); // Denormalize
        vec4 curValF = rval + value;    // Add
        curValF.rgb /= curValF.a;       // Renormalize
        newVal = convVec4ToRGBA8(curValF);
        ++numIterations;
以上步骤完成了模型的体素化。
 
五、通过mipmap来构建体素八叉树
 
当我们将模型体素化完毕后，我们需要使用一种数据结构对其进行管理，八叉树的主要思想是将一个3维空间（立方体）等分成八个次级立方体，并以此不停的迭代下去，直到每个立方体唯一代表一个物体（或其一部分）——叶子结点，所以其应该是一个从根节点进行扩展的树形结构，但是现在的情况是我们已经完成了体素化——即所有的叶子节点已经被计算出来，那么通过将相邻的八个体素合并成一个，并不断的迭代直到合并成一个和模型的AABB包围盒同样大小的体素来构建八叉树——是不是很像mipmap呢。
 
1、生成mipmap之前，先计算一下漫反射吧。
 
VCTRenderer的作者在生成mipmap之前先计算了漫反射和阴影，且分成了直接光和经过一次反射的间接光两部分的漫反射计算。
 
第一步计算体素的可见度（即阴影），这里作者采用了两种方式，第一种是根据计算得到的阴影贴图进行计算，其代码就是三、3：Exponential Variance Shadow Map中所述的代码。
 
第二种方式则是使用了光线步进——ray marching的方法来计算体素的可见度：
 
// 光线追踪计算阴影
float TraceShadow(vec3 position, vec3 direction, float maxTracingDistance) 
    // 缩放系数
    float k = traceShadowHit * traceShadowHit;
    // 计算体素在体素空间中的大小
    float voxelTexSize = 1.0f / volumeDimension;
    // 初始时移动两倍于体素大小的距离来避免发生自碰撞
    float dst = voxelTexSize * 2.0f;
    vec3 samplePos = direction * dst + position;
    // control variables
    float visibility = 0.0f;
    // accumulated sample
    float traceSample = 0.0f;
    while (visibility <= 1.0f && dst <= maxTracingDistance) 
	//判断边界条件
        if (samplePos.x < 0.0f || samplePos.y < 0.0f || samplePos.z < 0.0f || samplePos.x > 1.0f || samplePos.y > 1.0f || samplePos.z > 1.0f) 
            break; 
	//ceil返回大于等于X的最小整数值
        traceSample = ceil(texture(voxelAlbedo, samplePos).a) * k;
        // 如果透明度乘上缩放系数大于1则返回不可见
        if(traceSample > 1.0f - EPSILON) 
	    return 0.0f; 
        // 计算可见度
        visibility += (1.0f - visibility) * traceSample / dst;
        // 按光照方向移动
        dst += voxelTexSize;
        samplePos = direction * dst + position;
    return 1.0f - visibility;
第二步则是光照的计算，这里VCTRenderer的作者虽然函数名字写的是BRDF但是与我所知的BRDF相差甚远，更像是针对体素的立方体形状的一种trick算法（个人猜想，BRDF中的D项（法线分布函数）和G项（几何遮蔽函数）在这里并没有很大的作用，因为我们是按照体素进行光照计算，法线分布和几何遮蔽自然是根据体素的堆叠方式的不同而不同，）：
 
vec3 BRDF(Light light, vec3 normal, vec3 albedo)
    float nDotL = 0.0f;
    if(normalWeightedLambert == 1)
        vec3 weight = normal * normal;
        // 对X、Y、Z轴三个方向分别计算与光线的夹角
        float rDotL = dot(vec3(1.0, 0.0, 0.0), light.direction);
        float uDotL = dot(vec3(0.0, 1.0, 0.0), light.direction);
        float fDotL = dot(vec3(0.0, 0.0, 1.0), light.direction);
        rDotL = normal.x > 0.0 ? max(rDotL, 0.0) : max(-rDotL, 0.0);
        uDotL = normal.y > 0.0 ? max(uDotL, 0.0) : max(-uDotL, 0.0);
        fDotL = normal.z > 0.0 ? max(fDotL, 0.0) : max(-fDotL, 0.0);
        // 按照权重得到夹角大小
        nDotL = rDotL * weight.x + uDotL * weight.y + fDotL * weight.z;
        nDotL = max(dot(normal, light.direction), 0.0f);
    return light.diffuse * albedo * nDotL;
return vec4(BRDF(light, normal, albedo) * visibility, visibility);
 
return vec4(BRDF(light, normal, albedo) * falloff * visibility, visibility);
 
 return vec4(BRDF(light, normal, albedo) * falloff * spotFalloff * visibility, visibility);
 
2、间接光的漫反射计算（注意需要先生成Mipmap，我这里稍微调整了一下顺序）
 
这里使用了cone trace算法，每个圆锥的示意图如下：
 
由上图可知，圆锥体并不是真正意义上的圆锥体，而是由不同level的3D纹理（或者说体素）并接而成的类锥体的形状。
 
每个圆锥体都需要起始点C0、方向Cd、角度、追踪的长度t，由这四个变量可以计算出：
 
由d可以计算出需要在哪个level的mipmap上进行采样：
 
Vsize代表的含义是体素在最高的level的mipmap时的大小（这里一般是128）。
 
使用单个圆锥体并不能满足对整个半球空间进行积分运算（光照计算的基础，每个点都要计算其法线半球内所有光线对该点的影响才能得到最后的光照结果）这个条件，因此需要多个圆锥体进行计算，圆锥体的数量越多、角度越小结果越趋近于真正的半球积分效果，VCTRenderer的作者这里使用了四个圆锥体进行光照计算。
 
vec4 CalculateIndirectLighting(vec3 position, vec3 normal)
    // 沿着法线前进2/volumeDimension单位距离，因为0级的mipmap并不是256^3而是128^3
    position = position + normal * (1.0f / (volumeDimension / 2.0f));
    vec4 diffuseTrace = vec4(0.0f);
    // 设定上向量——用于计算TBN的辅助向量
    const float aperture = 1.0f;
    vec3 guide = vec3(0.0f, 1.0f, 0.0f);
    if (abs(dot(normal, guide)) == 1.0f)
        guide = vec3(0.0f, 0.0f, 1.0f);
    // 计算切线和副切线
    vec3 right = normalize(guide - dot(normal, guide) * normal);
    vec3 up = cross(right, normal);
    for(int i = 0; i < 4; i++)
        vec3 coneDirection = normal;
        coneDirection += propagationDirections[i].x * right + propagationDirections[i].z * up;    
        coneDirection = normalize(coneDirection);      //确定光线追踪的方向
        diffuseTrace += TraceCone(position, coneDirection, aperture) * diffuseConeWeights[i];
    return clamp(diffuseTrace, 0.0f, 1.0f);
vec4 TraceCone(vec3 position, vec3 direction, float aperture)
    uvec3 visibleFace;
    float anisoDimension = volumeDimension / 2.0f;
    // 只对看得见的面的mipmap进行光线追踪
    visibleFace.x = (direction.x < 0.0) ? 0 : 1;
    visibleFace.y = (direction.y < 0.0) ? 2 : 3;
    visibleFace.z = (direction.z < 0.0) ? 4 : 5;
    // 每个轴方向上的权重
    vec3 weight = direction * direction;
    float voxelSize = 1.0f / anisoDimension;
    // 为了避免自身碰撞移动一个体素的距离
    float dst = voxelSize;
    float diameter = aperture * dst;
    vec3 samplePos = position + direction * dst;
    float mipLevel = 0.0f;
    vec4 coneSample = vec4(0.0f);
    vec4 anisoSample = vec4(0.0f);
    if(samplePos.x < 0.0f || samplePos.y < 0.0f || samplePos.z < 0.0f || samplePos.x > 1.0f || samplePos.y > 1.0f || samplePos.z > 1.0f) 
        return coneSample; 
    //迭代光线追踪
    while(coneSample.a <= 1.0f && dst <= maxTracingDistanceGlobal)
        if (checkBoundaries > 0 && (samplePos.x < 0.0f || samplePos.y < 0.0f || samplePos.z < 0.0f || samplePos.x > 1.0f || samplePos.y > 1.0f || samplePos.z > 1.0f)) 
            break; 
        // 根据光线前进的距离计算需要采样的mip等级
        mipLevel = log2(diameter * anisoDimension);
        mipLevel = max(mipLevel - 1.0f, 0.0f);
        // 根据miplevel进行各向异性采样
        anisoSample = weight.x * textureLod(voxelTexMipmap[visibleFace.x], samplePos, mipLevel)
                    + weight.y * textureLod(voxelTexMipmap[visibleFace.y], samplePos, mipLevel)
                    + weight.z * textureLod(voxelTexMipmap[visibleFace.z], samplePos, mipLevel);
        coneSample += (1.0f - coneSample.a) * anisoSample;
        // 光线步进
        dst += max(diameter, voxelSize);
        diameter = dst * aperture;
        samplePos = direction * dst + position;
    return coneSample;
3、生成mipmap
 
VCTRenderer的作者将其分成了两个步骤，第一步是从原始的3D纹理中生成第0级的mipmap（并没有把原始的3D纹理作为最终的体素的mipmap纹理）。
 
第二步则是根据第一步生成的Level0的3Dmipmap纹理生成其余level的mipmap，shader十分简单：
 
#version 430
layout (local_size_x = 8, local_size_y = 8, local_size_z = 8) in;
layout(binding = 0, rgba8) uniform image3D voxelMipmapDst[6];
layout(binding = 5) uniform sampler3D voxelMipmapSrc[6];
uniform vec3 mipDimension;
uniform int mipLevel;
const ivec3 anisoOffsets[] = ivec3[8]
	ivec3(1, 1, 1),
	ivec3(1, 1, 0),
	ivec3(1, 0, 1),
	ivec3(1, 0, 0),
	ivec3(0, 1, 1),
	ivec3(0, 1, 0),
	ivec3(0, 0, 1),
	ivec3(0, 0, 0)
void FetchTexels(ivec3 pos, int dir, inout vec4 val[8]) 
	for(int i = 0; i < 8; i++)
		val[i] = texelFetch(voxelMipmapSrc[dir], pos + anisoOffsets[i], mipLevel);
void main()
	if(gl_GlobalInvocationID.x >= mipDimension.x || gl_GlobalInvocationID.y >= mipDimension.y || gl_GlobalInvocationID.z >= mipDimension.z) 
		return;
	ivec3 writePos = ivec3(gl_GlobalInvocationID);
	ivec3 sourcePos = writePos * 2;
	// fetch values
	vec4 values[8];
	// x -
	FetchTexels(sourcePos, 0, values);
	imageStore(voxelMipmapDst[0], writePos, 
		values[0] + values[4] * (1 - values[0].a) + 
		values[1] + values[5] * (1 - values[1].a) +
		values[2] + values[6] * (1 - values[2].a) +
		values[3] + values[7] * (1 - values[3].a)) * 0.25f
	// x +
	FetchTexels(sourcePos, 1, values);
    imageStore(voxelMipmapDst[1], writePos, 
		values[4] + values[0] * (1 - values[4].a) +
    	values[5] + values[1] * (1 - values[5].a) +
    	values[6] + values[2] * (1 - values[6].a) +
    	values[7] + values[3] * (1 - values[7].a)) * 0.25f
	// y -	
	FetchTexels(sourcePos, 2, values);
    imageStore(voxelMipmapDst[2], writePos, 
		values[0] + values[2] * (1 - values[0].a) +
    	values[1] + values[3] * (1 - values[1].a) +
    	values[5] + values[7] * (1 - values[5].a) +
    	values[4] + values[6] * (1 - values[4].a)) * 0.25f
	// y +
	FetchTexels(sourcePos, 3, values);
    imageStore(voxelMipmapDst[3], writePos, 
		values[2] + values[0] * (1 - values[2].a) +
    	values[3] + values[1] * (1 - values[3].a) +
    	values[7] + values[5] * (1 - values[7].a) +
    	values[6] + values[4] * (1 - values[6].a)) * 0.25f
	// z -
	FetchTexels(sourcePos, 4, values);
    imageStore(voxelMipmapDst[4], writePos, 
		values[0] + values[1] * (1 - values[0].a) +
    	values[2] + values[3] * (1 - values[2].a) +
    	values[4] + values[5] * (1 - values[4].a) +
    	values[6] + values[7] * (1 - values[6].a)) * 0.25f
	// z +
	FetchTexels(sourcePos, 5, values);
    imageStore(voxelMipmapDst[5], writePos, 
		values[1] + values[0] * (1 - values[1].a) +
    	values[3] + values[2] * (1 - values[3].a) +
    	values[5] + values[4] * (1 - values[5].a) +
    	values[7] + values[6] * (1 - values[7].a)) * 0.25f
六、全局光照
 
这里的全局光照计算其算法与五中计算直接漫反射和间接漫反射的算法差别不大，这里主要说一下几个没叙述过的点。
 
1、通过深度和纹理坐标计算出世界坐标
 
因为纹理坐标覆盖的是整个窗口，所以本质上此时纹理坐标等价于顶点转换到ndc空间后*2-1，因此可由下述函数计算出世界坐标：
 
// 所以texCoord可以变换到[-1，1]范围内来当作ndc空间内该片段所对应的点的x、y坐标
// 随后通过采样深度贴图获得深度值，再乘上逆矩阵便可以得到坐标
vec3 PositionFromDepth()
    float z = texture(gDepth, texCoord).x * 2.0f - 1.0f;
    vec4 projected = vec4(texCoord * 2.0f - 1.0f, z, 1.0f);
    projected = inverseProjectionView * projected;
    return projected.xyz / projected.w;
2、射线与AABB求交
 
Origin代表射线的起始点，dir是射线方向。
 
normal是平面的法线，X是平面上一点，d是平面到原点的距离。
 
由此求出最大的t并保证交点在AABB上，即可判定射线与AABB的哪个面相交。
 
bool IntersectRayWithWorldAABB(vec3 ro, vec3 rd, out float enter, out float leave)
    //因为AABB的法线都是(0,1,0),(1,0,0),(0,0,1)
    //所以直接对应分量相除即可
    vec3 tempMin = (worldMinPoint - ro) / rd; 
    vec3 tempMax = (worldMaxPoint - ro) / rd;
    vec3 v3Max = max (tempMax, tempMin);
    vec3 v3Min = min (tempMax, tempMin);
    leave = min (v3Max.x, min (v3Max.y, v3Max.z));
    enter = max (max (v3Min.x, 0.0), max (v3Min.y, v3Min.z));    
    return leave > enter;
3、直接光光照计算
 
vec3 BRDF(Light light, vec3 N, vec3 X, vec3 ka, vec4 ks)
    vec3 L = light.direction;
    vec3 V = normalize(cameraPosition - X);
    vec3 H = normalize(V + L);
    float dotNL = max(dot(N, L), 0.0f);
    float dotNH = max(dot(N, H), 0.0f);
    float dotLH = max(dot(L, H), 0.0f);
    // 解码高光功率
    float spec = exp2(11.0f * ks.a + 1.0f);
    // 菲涅尔效应
    vec3 fresnel = ks.rgb + (1.0f - ks.rgb) * pow(1.0f - dotLH, 5.0f);
    // 高光因子
    float blinnPhong = pow(dotNH, spec);
    // 能量守恒
    blinnPhong *= spec * 0.0397f + 0.3183f;
    // 高光
    vec3 specular = ks.rgb * light.specular * blinnPhong * fresnel;
    // 漫反射
    vec3 diffuse = ka.rgb * light.diffuse;
    return (diffuse + specular) * dotNL;
return BRDF(light, normal, position, albedo, specular) * visibility;
 
return BRDF(light, normal, position, albedo, specular) * falloff * visibility;
 
return BRDF(light, normal, position, albedo, specular) * falloff * spotFalloff * visibility;
 
4、间接光光照计算
 
vec4 CalculateIndirectLighting(vec3 position, vec3 normal, vec3 albedo, vec4 specular, bool ambientOcclusion)
    vec4 specularTrace = vec4(0.0f);
    vec4 diffuseTrace = vec4(0.0f);
    vec3 coneDirection = vec3(0.0f);
    if(any(greaterThan(specular.rgb, specularTrace.rgb)))
        vec3 viewDirection = normalize(cameraPosition - position);
        vec3 coneDirection = reflect(-viewDirection, normal);
        coneDirection = normalize(coneDirection);
        // 镜面反射锥体角度大小设置，最低数值为1单位梯度，再小的话会严重影响性能
        float aperture = clamp(tan(HALF_PI * (1.0f - specular.a)), 0.0174533f, PI);
        specularTrace = TraceCone(position, normal, coneDirection, aperture, false);
        specularTrace.rgb *= specular.rgb;
    if(any(greaterThan(albedo, diffuseTrace.rgb)))
        // 这里的漫反射调整了角度的大小以及圆锥体的数量
        const float aperture = 0.57735f;
        vec3 guide = vec3(0.0f, 1.0f, 0.0f);
        if (abs(dot(normal,guide)) == 1.0f)
            guide = vec3(0.0f, 0.0f, 1.0f);
        // Find a tangent and a bitangent
        vec3 right = normalize(guide - dot(normal, guide) * normal);
        vec3 up = cross(right, normal);
        for(int i = 0; i < 6; i++)
            coneDirection = normal;
            coneDirection += diffuseConeDirections[i].x * right + diffuseConeDirections[i].z * up;
            coneDirection = normalize(coneDirection);
            // cumulative result
            diffuseTrace += TraceCone(position, normal, coneDirection, aperture, ambientOcclusion) * diffuseConeWeights[i];
        diffuseTrace.rgb *= albedo;
    vec3 result = bounceStrength * (diffuseTrace.rgb + specularTrace.rgb);
    return vec4(result, ambientOcclusion ? clamp(1.0f - diffuseTrace.a + aoAlpha, 0.0f, 1.0f) : 1.0f);
以上就是对VCTRenderer的大体解析，对于光照部分一些细节处理，我并没有完全看懂，作者是按照Interactive Indirect Illumination Using Voxel Cone Tracing这篇论文写得cone trace的代码，如果有疑问可以去看一下这篇论文。
				
这是一个用于研究现代图形技术的简单DX11.1渲染引擎。 特别是根据Crassin等人的“使用Voxel圆锥体跟踪进行交互式间接照明”的Voxel圆锥体跟踪实现。 （Cyril Crassin，Fabrice Neyret，Miguel Saintz，Simon Green和Elmar Eisemann） 
 易于使用。 它具有基本的obj / mtl加载器，支持静态场景，使用effect11框架和预编译的fx着色器。 我添加了抽象注释并绘制了几种方案来简化代码读取。 如果您想编程任何技术，可以从DefaultShader和RenderTick开始。
 要求：Microsoft Redistributable 2013，d3dcompiler_47.dll，DX11.1兼容适配器（或至少DX10.0兼容适配器来运行应用程序）。
 演示：  : 说明：  : 
 作者：Dontsov Val
				
之前很早就看到了UE4中的基于Sparse Voxel Octree的RTGI，效果很酷，一直尝试作些研究与实现，但苦于没机会。前段得空，抽时间学习了一下，这里小结一下备忘。
整个算法主要分类几个过程：体素化、Mipmap OCTree、Cone Tracing。
1. Voxelization
体素化整个GI算法的基础。这里体素化可以采用的方法也比较多，主要有以下几种：
				
废话就不多说了，开始。。。
    之前很早就看到了UE4中的基于Sparse Voxel Octree的RTGI，效果很酷，始终尝试作些研讨与实现，但苦于没机会。前段无暇，抽时间学习了一下，这里小结一下备忘。
    整个算法重要分类几个进程：体素化、Mipmap OCTree、Cone Tracing。
1. Voxelization
    体素化整个GI算法的基本。这里...
					
体素锥形追踪（Voxel Cone Tracing）在实时全局光照中是如何利用DirectCompute和OpenGL技术提高渲染效率的？请结合具体技术细节进行说明。
				
本篇文章主要是实现python 自然语言处理包 gensim 中用于词向量建模的 word2vec算法。示例代码如下：# encoding=utf-8import logging
import sysfrom gensim.models import Word2Vecif __name__ == '__main__':
    logging.basicConfig(format='%(asctim
				
本篇文章主要是实现Python 自然语言处理包 gensim 中用于长文本向量建模的  doc2vec算法。示例代码如下：#!/usr/bin/env python3
# -*- coding: utf-8 -*-import logging
import multiprocessing
import os.path
import sysfrom gensim import utils
from g