# How to implement PCSS properly?

I'm trying to implement PCSS in OpenGL/GLSL, but I have problems understanding the details especially the conversions between the coordinate spaces. The existing implementations are so different from each other and I'm just more confused. I'm sure they all good implementations I just don't understand them.

## My implementation

In the vertex shader I don't do anything special, I just compute positions in different coordinate spaces.

#version 460 core

layout(location=0) in vec3 i_position;
layout(location=1) in vec3 i_normal;

uniform mat4 u_model;
uniform mat4 u_view;
uniform mat4 u_projection;

uniform mat4 u_light_view;
uniform mat4 u_light_projection;

out vec3 io_position;
out vec3 io_normal;
out vec4 io_lvs_position;
out vec4 io_lcs_position;

void main(){
io_position = vec3(u_model * vec4(i_position, 1.0));            //world space position
gl_Position = u_projection * u_view * vec4(io_position, 1.0);   //clip space position
io_normal = mat3(u_model) * i_normal;                           //world space normal
io_lvs_position = u_light_view * vec4(io_position, 1.0);        //light view space position
io_lcs_position = u_light_projection * io_lvs_position;         //light clip space position
}


float compute_pcss_shadow(){
float search_region_width = compute_search_region_width();
float average_blocker_depth = compute_average_blocker_depth(search_region_width);
if(average_blocker_depth == -1.0){
return 1.0;
}
float blocker_distance = compute_blocker_distance(average_blocker_depth);
float penumbra_width = compute_penumbra_width(blocker_distance);
float pcf_width = compute_pcf_width(penumbra_width);
return compute_pcss(pcf_width);
}


I imagine the compute_search_region_width function something like this:

We basically project the light's width into the near plane and we assume everything inside the triangle is a blocker. Because of the 2 similar triangles search_region_width / u_light_width and (lvs_distance - u_near_plane) / lvs_distance are equal, so I can compute search_region_width. However 3 out of 4 implementations (except diharaw's) divide the light's width by the frustum's width. I think they do that to convert the width from light view space to the (0;1) interval because in the next function we'll use this value in texture space. But if it's true it makes more sense to me computing the width in light view space and then dividing, rather than dividing the light's width. This is my implementation:

float compute_search_region_width(){
float lvs_distance = -io_lvs_position.z;  //shaded point's distance from the light in light view space
return u_light_size * (lvs_distance - u_near_plane) / lvs_distance / u_frustum_width;
}


Am I right or wrong with that u_frustum_width division?

The compute_average_blocker_depth function is just sampling the shadow map in the search region and averaging the blockers' depth.

float compute_average_blocker_depth(float search_region_width){
vec3 uv = io_lcs_position.xyz / io_lcs_position.w;  //conversion from clip space to NDC
uv = uv * 0.5 + 0.5;                                //conversion from the (-1;1) interval to the (0;1)
float real_depth = uv.z;
int blocker_count = 0;
float blocker_depth_sum = 0;

int step_count = u_kernel_size / 2;
float step_uv = search_region_width / step_count;
for(float x=-step_count;x<=step_count;x++){
for(float y=-step_count;y<=step_count;y++){
vec2 offset = vec2(x, y) * step_uv;
float depth = texture(u_shadow_map, uv.xy + offset).r;
if(depth < real_depth){
blocker_count++;
blocker_depth_sum += depth;
}
}
}
if(blocker_count == 0){
return -1.0;
}else{
return blocker_depth_sum / blocker_count;
}
}


In the compute_blocker_distance function we convert the average blocker depth from the (0;1) interval to light view space. My implementation is pretty simple and similar to diharaw's and asylum2010's implementation.

float compute_blocker_distance(float average_blocker_depth){
return average_blocker_depth * (u_far_plane - u_near_plane) + u_near_plane;
}


However in pboechat's implementation this step is completely missing and in the Nvidia implementation it looks like this:

float ZClipToZEye(float zClip)
{
return g_LightZFar*g_LightZNear / (g_LightZFar - zClip*(g_LightZFar-g_LightZNear));
}


This piece of code doesn't make any sense to me, maybe this is because this is the only HLSL code.

In the compute_penumbra_width function we determine the shadow penumbra's width in light view space, again based on similar triangles. This is how I imagine it:

Because penumbra_width / u_light_width and (lvs_distance - average_blocker_depth) / average_blocker_depth are equal, we can compute penumbra_width in light view space. My implementation is similar to asylum2010' and the Nvidia implementation:

float compute_penumbra_width(float blocker_distance){
float lvs_distance = -io_lvs_position.z; //shaded point's distance from the light in light view space
return u_light_size * (lvs_distance - blocker_distance) / blocker_distance / u_frustum_width;
}


I'm not entirely sure why we divide by the frustum width, but later we'll use the result in texture space rather than in light view space, so maybe this is that conversion again. In pboechat's implementation we use the shaded point's depth in the (0;1) interval, not the distance in light view space, which doesn't make any sense to me. In diharaw's implementation the computation is something like this: abs(zReceiver - zBlocker) / zBlocker and I don't understand it.

The compute_pcf_width is the function what I don't understand at all. At this point I think penumbra_width is in texture space because of the division by u_frustum_width. So what's the point of more conversions? But this is the point where the other implementations are very similar. Multiply the penumbra's width with the near plane and then divide by the shaded point's light view space distance (in pboechat's implementation divide by the depth). My implementation is like this:

float compute_pcf_width(float penumbra_width){
vec3 uv = io_lcs_position.xyz / io_lcs_position.w;  //conversion from clip space to NDC
uv = uv * 0.5 + 0.5;                                //conversion from the (-1;1) interval to the (0;1)
float real_depth = uv.z;
return penumbra_width * u_near_plane / real_depth;
}


I'm using the depth and not the light view space distance, because it looks better, but I don't understand the multiplication and the division in this function. So what is the point of the compute_pcf_width function and how to implement it properly?

The last function, compute_pcss is a variable kernel sized PCF filter.

float compute_pcss(float pcf_width){
vec3 uv = io_lcs_position.xyz / io_lcs_position.w;
uv = uv * 0.5 + 0.5;
float real_depth = uv.z;
float result = 0;

int step_count = u_kernel_size / 2;
float step_uv = pcf_width / step_count;
for(float x=-step_count;x<=step_count;x++){
for(float y=-step_count;y<=step_count;y++){
vec2 offset = vec2(x, y) * step_uv;
float depth = texture(u_shadow_map, uv.xy + offset).r + bias;
if(depth > real_depth){
result += 1.0;
}else{
result += 0.5;
}
}
}
return result / (u_kernel_size * u_kernel_size);
}


This is the result with u_kernel_size = 27:

The shadow is getting blurrier with distance, but the edge of the shadow looks a little strange. Also you can see some sort of issue on the right side, maybe blurred shadow acne. So if you read all this, please tell me what is the correct implementation of PCSS and why.