# Screen space AO technique that is faster than HBAO and just as accurate

doc.md

Planar Screen Space Ambient Occlusion

PSSAO works off of SSAO, but also includes how much a depth difference matters in the end result. Let's say a is our depth/distance value. This distance is: (in pseudocode) a = (dot(centerPixelNormal, neighboringPixel3DPos) - dot(centerPixelNormal, centerPixel3DPos)) / distance( neighboringPixelPos * cameraDimensions, centerPixelPos * cameraDimensions)

Now put it in this formula:

1/(1+a)

Where a is the aforementioned depth difference. The result is the occlusion from that pixel. Now, the occlusion from an array of pixels in a direction is the maximum occlusion value. And then sample a circle kernel for the ambient occlusion, and sum the occlusion from every direction. The result is to be put in the ambient occlusion buffer.

Why This Works

The First Part: Why the Max Occlusion Value Works

You may be wondering, why I say when calculating the occlusion in a direction, to calculate the maximum occlusion value in that direction to accumulate the occlusion value. This is because of perspective.

So?

Say in one direction, you have two occluding pillars A and B.

This image is what the pillars look like, looking from the perspective of the center of the hemisphere.

A is taller than B, but is in the same direction relative to the circle kernel cast out by the algorithm (See in my explanation of the algo, in the end.) So let's say A occludes 50% of that direction. B occludes 20%. If B is removed or behind A, the resulting occlusion is 50% (Because A is the only pillar in that direction.)

If we use summation to calculate the occlusion, we get 70%. But if we look at the original image, and the image where A is the only pillar visible, we see that the amount of occlusion is the same. And since A is more occluding than B, that means that the max occlusion of all the occlusion percentages is the amount of occlusion.

Second Part: Why we divide by distance from the center pixel without including depth

Say you have a pillar A that occludes 50%. If you move that pillar far away, will it occlude the hemisphere? Most likely not.

GLSL Code:

uniform vec2 cameraDimensions;
uniform vec2 cameraPixelSize;
uniform mat4 projectionMatrix;
uniform mat4 viewMatrix;
uniform sampler2D depthTex;

// https://sibaku.github.io/computer-graphics/2017/01/10/Camera-Ray-Generation.html
vec3 getPos(ivec2 px, float depth, mat4 PInv, mat4 VInv)
{
    vec2 pxNDS = (vec2(px) / cameraPixelSize) * 2. - 1.;
    // z = -1 equals a point on the near plane, i.e. the screen
    vec3 pointNDS = vec3(pxNDS, -1.);
    vec4 pointNDSH = vec4(pointNDS, 1.0);
    vec4 dirEye = PInv * pointNDSH;
    dirEye.w = 0.;
    vec3 dirWorld = (VInv * dirEye).xyz;
    return normalize(dirWorld) * depth;
}

// https://www.shadertoy.com/view/fsVczR
vec3 computeNormalImproved( sampler2D depth, ivec2 p, mat4 PInv, mat4 VInv )
{
    float c0 = texelFetch(depth,p           ,0).r;
    float l2 = texelFetch(depth,p-ivec2(2,0),0).r;
    float l1 = texelFetch(depth,p-ivec2(1,0),0).r;
    float r1 = texelFetch(depth,p+ivec2(1,0),0).r;
    float r2 = texelFetch(depth,p+ivec2(2,0),0).r;
    float b2 = texelFetch(depth,p-ivec2(0,2),0).r;
    float b1 = texelFetch(depth,p-ivec2(0,1),0).r;
    float t1 = texelFetch(depth,p+ivec2(0,1),0).r;
    float t2 = texelFetch(depth,p+ivec2(0,2),0).r;
    
    float dl = abs(l1*l2/(2.0*l2-l1)-c0);
    float dr = abs(r1*r2/(2.0*r2-r1)-c0);
    float db = abs(b1*b2/(2.0*b2-b1)-c0);
    float dt = abs(t1*t2/(2.0*t2-t1)-c0);
    
    vec3 ce = getPos(p,c0,PInv,VInv);

    vec3 dpdx = (dl<dr) ?  ce-getPos(p-ivec2(1,0),l1,PInv,VInv) : 
                          -ce+getPos(p+ivec2(1,0),r1,PInv,VInv) ;
    vec3 dpdy = (db<dt) ?  ce-getPos(p-ivec2(0,1),b1,PInv,VInv) : 
                          -ce+getPos(p+ivec2(0,1),t1,PInv,VInv) ;

    return normalize(cross(dpdx,dpdy));
}

float calcOcclusionFromPixel(vec2 pos, vec2 centerPos, float centerPixDepth, vec3 centerPixNormal, mat4 PInv, mat4 VInv)
{
    float neighborPixelDepth = texture(depthTex, pos).r;
    vec3 centerPixel3DPos = getPos(ivec2(centerPos * cameraPixelSize), centerPixDepth, PInv, VInv);
    vec3 neighborPixel3DPos = getPos(ivec2(pos * cameraPixelSize), neighborPixelDepth, PInv, VInv);
    float a = dot(centerPixNormal, neighborPixel3DPos) - dot(centerPixNormal, centerPixel3DPos);
    if(a <= 0.0) return 0.0;
    a /= 1.0 + distance(centerPixel3DPos, neighborPixel3DPos);
    a = 1.0 / (1.0 + a);
    return 1.0 - a;
}

float calcOcclusionInDirection(int steps, vec2 direction, vec2 centerPos, float centerPixDepth, vec3 centerPixNormal, mat4 PInv, mat4 VInv, out float toInclude)
{
    float max_ = -100000.0;
    for(int i = 0; i != steps; i++)
    {
        vec2 pix = floor(direction * float(i)) / cameraPixelSize;
        pix += centerPos;
        if( all(equal(pix, clamp(pix, 0.0, 1.0))) ) toInclude += 1.0;
        else continue;
        max_ = max(max_, calcOcclusionFromPixel(pix, centerPos, centerPixDepth, centerPixNormal, PInv, VInv));
    };
    
    return max_;
}

uniform int radiusSteps;
uniform int angularSteps;

float calcOcclusion(vec2 centerPos)
{
    float sum = 0.0;
    float toInclude = 0.0;
    float centerPixDepth = texture(depthTex, centerPos).r;
    mat4 PInv = inverse(projectionMatrix);
    mat4 VInv = inverse(viewMatrix);
    vec3 centerPixNormal = computeNormalImproved(depthTex, ivec2(centerPos * cameraPixelSize), PInv, VInv);
    
    for(int i = 0; i != angularSteps; i++)
    {
        float angle = (float(i) / float(angularSteps)) * 3.141592653 * 2.0;
        sum += calcOcclusionInDirection(radiusSteps, normalize(vec2(cos(angle), sin(angle))), centerPos, centerPixDepth, centerPixNormal, PInv, VInv, toInclude);
    };
    
    return sum / toInclude;
}