RealTime-Rendering15-GPU Culling

发表于 更新于 分类于 ， ， 阅读次数： 本文字数： 1.9k 阅读时长 ≈ 7 分钟

GPU Culling

基于computeShader实现GPU的视锥剔除和遮挡剔除.

frustumCulling

GPU提出的第一步先判断包围盒是否完全位于视锥体外部，需要构造frustum的六个平面以及场景包围盒，循环frustum的6个面分别与包围盒求交，如果包围盒在任意一个平面的外部，则直接剔除。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
bool frustumCulling(mat4 mvp, BoundingBox bbox)
{
  vec3 center = (bbox.min.xyz + bbox.max.xyz) * 0,5;
  vec3 extents = bbox.max.xyz - center;
  for(int i = 0; i < 6; ++i)
  {
    vec3 normal = frustumPlanes[i].xyz;
    float d = frustumPlanes[i].w;
    float radius = dot(extents, abs(normal));
    float distance = dot(center, normal) + d;
    if(distance < -radius)
    {
      return true;
    }
  }
  return false;
}

extract Planes From Matrix

而要获取构成frustum的六个面, 首先需要从viewprojection矩阵中提取这六个面.假设世界空间的某个向量$\vec{v} = (x,y,z,w)$, 矩阵$M = m_{ij}$是viewProjectionMatrix，使用$M$变换$\vec{v}$可以表示为矩阵和向量的乘法：

01.png

其中 • 表示向量的点乘，$Row_i=(m_{i1},m_{i2}, m_{i3}, m_{i4})$表示矩阵$M$的第$i$行。

经过此变换后，顶点$v^′$位于齐次裁剪空间中。在此空间中，视锥体实际上是一个轴对齐的盒状区域，其尺寸由具体API决定。若顶点$v^′$位于该盒内，则未经变换的顶点$v$便处于”未经变换的”视锥体内部。在OpenGL中，当$v^′$的分量满足以下不等式时，$v^′$位于该盒内：

$−w^′<x^′<w^′$

$−w^′<y^′<w^′$

$−w^′<z^′<w^′$

条件 结论
−w′<x′ x′位于左裁剪平面的内侧半空间
x′<w′ x′位于右裁剪平面的内侧半空间
−w′< y′ y′位于下裁剪平面的内侧半空间
y′<w′ y′位于上裁剪平面的内侧半空间
−w′<z′ z′位于近裁剪平面的内侧半空间
z′<w′ z′位于远裁剪平面的内侧半空间

现在假设我们要测试$x^′$是否位于左裁剪平面的内侧半空间。需要满足以下不等式：

$−w^′<x^′$

利用本章开头的信息，我们可以将该不等式重写为：

$−(\vec{v}⋅Row_4)<(\vec{v}⋅Row_1)$

这等价于：

$0<(\vec{v}⋅Row_4)+(\vec{v}⋅Row_1)$

最终得到：

$0<\vec{v}⋅(Row_4+Row_1)$

这实际上已经代表了“未经变换的”视锥体左裁剪平面的平面方程：

02.png

w等于1，平面方程可以进一步简化为：

03.png

这样就得到了标准平面方程：

04.png

其中：

05.png

至此，我们证明了视锥体的左裁剪平面可以直接从投影矩阵中提取。需要注意的是，所得平面方程是未归一化的（即平面的法向量不是单位向量），并且该法向量指向内侧半空间。这意味着，当顶点$v$位于左裁剪平面的内侧半空间时，有 $0<ax+by+cz+d$。

归一化平面的完整步骤:

  1. 计算法线长度

06.png

  1. 将平面方程的$abcd$都除以$|n|$

07.png

  1. 得到归一化后的平面方程：

08.png

code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
Frustum extractFrustum(const glm::mat4& viewProj) 
{
    Frustum frustum;
    glm::mat4 m = glm::transpose(viewProj); // 转为列访问
    
    std::array<glm::vec4, 6> coeffs = {
        m[3] + m[0],  // Left
        m[3] - m[0],  // Right
        m[3] + m[1],  // Bottom
        m[3] - m[1],  // Top
        m[3] + m[2],  // Near
        m[3] - m[2]   // Far
    };
    
    for (int i = 0; i < 6; ++i) 
    {
        float len = glm::length(glm::vec3(coeffs[i]));
        frustum.planes[i].normal = glm::vec3(coeffs[i]) / len;
        frustum.planes[i].distance = coeffs[i].w / len;
    }
    
    return frustum;
}

References:

ExtractPlanes

AABB & Plane intersectionTest

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
int TestAABBPlane(AABB b, Plane p) 
{
    // Convert AABB to center-extents representation
    Point c = (b.max + b.min) * 0.5f; // Compute AABB center
    Point e = b.max - c; // Compute positive extents

    // Compute the projection interval radius of b onto L(t) = b.c + t * p.n
    float r = e[0]*Abs(p.n[0]) + e[1]*Abs(p.n[1]) + e[2]*Abs(p.n[2]);

    // Compute distance of box center from plane
    float s = Dot(p.n, c) - p.d;

    // Intersection occurs when distance s falls within [-r,+r] interval
    return Abs(s) <= r;
}

aabb & plane intersection

OcclusionCulling

这里使用的遮挡剔除技术通过采样深度图比较判断深度的方式判断遮挡性，而场景中不同深度的物体包围盒对应的纹素空间比例是不一样的，离相机较远的物体一个像素可能会对应深度图的一片区域，因此，为了保证剔除的正确性，采用保守剔除策略，宁肯少剔除也不要错误剔除，需要事先对zBuffer手动生成mipmap，生成Hierarchy zBuffer，简写为HiZ。为什么不能使用硬件自动生成mipmap呢？因为硬件mipmap采用boxFilter的均值策略在此并不适用。这里为了防止误剔除，采用保守策略，mipmap要取最大值。

Generate HiZ

使用compute Shader实现手动mipmap：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
#version 460 core

layout(local_size_x = 16, local_size_y = 16) in;

// 输入纹理（上一级 Mip）
layout(binding = 0, r32f) uniform readonly image2D srcDepth;
// 输出纹理（当前级 Mip）
layout(binding = 1, r32f) uniform writeonly image2D dstDepth;

uniform ivec2 srcSize;
uniform ivec2 dstSize;

void main() 
{
    ivec2 dstCoord = ivec2(gl_GlobalInvocationID.xy);
    
    if (dstCoord.x >= dstSize.x || dstCoord.y >= dstSize.y) 
    {
        return;
    }
    
    // 计算源坐标（2x2 降采样）
    ivec2 srcCoord = dstCoord * 2;
    
    // 采样 2x2 区域，取最小深度（用于遮挡剔除，保留最近物体）
    // 如果是用于反射/折射，可能需要取最大深度
    float d0 = imageLoad(srcDepth, srcCoord).r;
    float d1 = imageLoad(srcDepth, srcCoord + ivec2(1, 0)).r;
    float d2 = imageLoad(srcDepth, srcCoord + ivec2(0, 1)).r;
    float d3 = imageLoad(srcDepth, srcCoord + ivec2(1, 1)).r;
    
    float maxDepth = max(max(d0, d1), max(d2, d3));
    
    imageStore(dstDepth, dstCoord, vec4(maxDepth, 0.0, 0.0, 0.0));
}

OcclusionCulling

有了HiZ就可以执行遮挡剔除了，首先第一步将AABB投影到屏幕空间

projectAABBToScreen

基于AABB的MinMax推导其8个顶点，将8个顶点依次投影到屏幕空间，计算xy的最小值最大值，即screenMin和screenMax，计算这个值的目的是为了推导aabb在纹理空间（屏幕空间）的覆盖范围，由此计算当前适合采样的HiZ层级.另外计算当前AABB的最小深度，用于后续深度判断：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
bool projectAABBToScreen(vec3 aabbMin, vec3 aabbMax, out vec2 screenMin, out vec2 screenMax, out float depthMin)
{
  vec3 corners[8];
  corners[0] = vec3(aabbMin.x, aabbMin.y, aabbMin.z);
  corners[1] = vec3(aabbMax.x, aabbMin.y, aabbMin.z);
  corners[2] = vec3(aabbMin.x, aabbMax.y, aabbMin.z);
  corners[3] = vec3(aabbMax.x, aabbMax.y, aabbMin.z);
  corners[4] = vec3(aabbMin.x, aabbMin.y, aabbMax.z);
  corners[5] = vec3(aabbMax.x, aabbMin.y, aabbMax.z);
  corners[6] = vec3(aabbMin.x, aabbMax.y, aabbMax.z);
  corners[7] = vec3(aabbMax.x, aabbMax.y, aabbMax.z);

  vec4 clipCorners[8];
  bool bAllBehind = true;
  screenMin = vec2(1.0, 1.0);
  screenMax = vec2(0.0, 0.0);
  depthMin = 1.0;

  for(int i = 0; i < 8; ++i)
  {
    clipCorners[i] = uViewProjectionMat * vec4(corners[i], 1.0);
    if(clipCorners[i].w > 0.0)
    {
      bAllBehind = false;
      vec3 ndc = clipCorners[i].xyz / clipCorners[i].w;
      vec2 screenPos = ndc.xyz * 0.5 + 0.5;
      screenMin = min(screenMin, screenPos);
      screenMax = max(screenMax, screenPos);
      depthMin = min(depthMin, ndc.z * 0.5 + 0.5);
    }
  }
  if(bAllBehind)
  {
    return false;
  }
  screenMin = clamp(screenMin, 0.0, 1.0);
  screenMax = clamp(screenMax, 0.0, 1.0);
  return true;
}

Calculate HiZ Level

上一步得到的screenMin和ScreenMax范围是[0, 1]，要想计算深度图的覆盖范围，还需要乘以纹理大小，然后取log2即可（参考mipmap）:

1
2
3
4
5
6
7
8
9
10
11
12
13
//AABB所有顶点都在相机后，直接剔除
if(projectAABBToScreen(aabbMin, aabbMax, screenMin, screenMax, depthMin))
{
  return true;
}
vec2 pixelSize = (screenMax - screenMin) * uViewportSize;
//小于1个像素的物体视为不可见，直接剔除
if(pixelSize.x < 1.0 && pixelSize.y < 1.0)
{
  return true;
}
float mipmapLevel = log2(max(pixelSize.x, pixelSize.y));
int mip = clamp(int(mipmapLevel), 0, uHizMipLevels - 1);

HiZ Occlusion Check

确定了HiZ层级，就可以对其进行深度采样了，采样点使用AABB投影到屏幕空间的矩形的四个顶点及中心点共五次采样：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
vec2 samples[5] = 
{
  vec2(screenMin.x, screenMin.y),
  vec2(screenMax.x, screenMin.y),
  vec2(screenMin.x, screenMax.y),
  vec2(screenMax.x, screenMax.y),
  vec2((screenMin.x + screenMax.x) * 0.5, (screenMin.y + screenMax.y) * 0.5)
};
float maxDepth = 0.0;
for(int i = 0; i < 5; ++i)
{
  float depth = textureLod(uHiZMap, samples[i], mip).r;
  //由于HiZMap是上一帧的数据，先将其反投影的世界空间，然后投影到当前帧，矫正深度值
  vec4 prevClipPos = vec4(samples[i] * 2.0 - 1.0, depth * 2.0 - 1.0, 1.0);
  vec4 worldPos = uPrevInvViewProj * prevClipPos;
  worldPos.xyz /= worldPos.w;

  vec4 currentClipPos = uViewProj * vec4(worldPos.xyz, 1.0);
  currentClipPos.xyz /= currentClipPos.w;
  float reprojectedDepth = currentClipPos.z + 0.1 / 9999.0f;
  maxDepth = max(maxDepth, reprojectedDepth);
  if(depthMin <= maxDepth)
  {
    return false;
  }
}
return true;

Generate IndirectDrawCommand

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
void main()
{
  int idx = gl_GlobalInvocationID.x;
  if (idx >= bboxes.length()) 
  {
    return;
  }
  BoundingBox bbox = bboxs[idx];
  if(frustumCulling(bbox))
  {
    return;
  }
  if(hizCulling(bbox))
  {
    return;
  }
  uint visibleIndex = atomicAdd(visibleCount, 1);
  outputCommands[visibleIndex] = inputCommands[idx];
  outputCommands[visibleIndex].baseInstance = idx;// 保持baseInstance不变，用于索引模型矩阵 
}