[path_tracing] implemented missing parts, fixed a lot of bugs, almost there

Author: PanosK92

data/shaders/restir_pt.hlsl

@@ -25,6 +25,7 @@ CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
 //==============================
 
 static const uint  INITIAL_CANDIDATE_SAMPLES   = 16;
+static const uint  LIGHT_RIS_CANDIDATE_SAMPLES = 4;
 static const float MIN_COS_AT_PRIMARY          = 1e-3f;
 static const float RUSSIAN_ROULETTE_PROB       = 0.85f;
 static const uint  RUSSIAN_ROULETTE_START      = 2;
@@ -338,12 +339,17 @@ float3 sample_sky(float3 dir)
 // direction is src's primary, so the stored radiance is view-dependent at rc w.r.t. the
 // reconnection shift, but the rc roughness gate bounds the resulting bias during reuse
 // max_bounces is the remaining path budget after this vertex (including this vertex's nee)
+// out_first_dir receives the first outgoing direction sampled at `start` so the caller can
+// populate rc_outgoing_dir for future variable-rc shift work; it is zeroed when no bounce
+// actually leaves start (max_bounces < 2 or first sample was rejected)
 float3 accumulate_subpath_radiance(
     PathPayload start,
     float3 start_view_dir,
     uint max_bounces,
-    inout uint seed)
+    inout uint seed,
+    out float3 out_first_dir)
 {
+    out_first_dir     = float3(0, 0, 0);
     float3 total      = float3(0, 0, 0);
     float3 throughput = float3(1, 1, 1);
 
@@ -360,6 +366,7 @@ float3 accumulate_subpath_radiance(
     float3      view_dir = start_view_dir;
     float       prev_brdf_pdf = 0.0f;
     float3      prev_normal   = start.hit_normal;
+    bool        first_dir_set = false;
 
     for (uint bounce = 1; bounce < max_bounces; bounce++)
     {
@@ -378,6 +385,13 @@ float3 accumulate_subpath_radiance(
         if (pdf < RESTIR_MIN_PDF || dot(nd, cur.hit_normal) <= 0.0f || any(isnan(nd)))
             break;
 
+        // capture the very first outgoing direction at the rc vertex
+        if (!first_dir_set)
+        {
+            out_first_dir = nd;
+            first_dir_set = true;
+        }
+
         float  unused_pdf;
         float3 brdf = evaluate_brdf(cur.albedo, cur.roughness, cur.metallic, cur.hit_normal, view_dir, nd, unused_pdf);
         throughput    *= brdf / pdf;
@@ -414,6 +428,106 @@ float3 accumulate_subpath_radiance(
     return total;
 }
 
+// builds a path sample candidate by directly sampling an analytical light or the sun cone
+// the candidate's rc vertex is the sampled light point (or sky direction for the sun) and
+// rc_radiance carries the emitted radiance toward the primary, the source_pdf returned is
+// in solid-angle measure at the primary so it can be combined with brdf-sampled candidates
+// in a single mixture-ris pool. point/spot delta lights are skipped: their solid-angle pdf
+// is effectively infinite and they are already handled by the spatial-pass direct lighting
+PathSample sample_light_candidate(
+    float3 primary_pos,
+    float3 primary_normal,
+    inout uint seed,
+    out float source_pdf)
+{
+    PathSample s;
+    s.rc_pos          = float3(0, 0, 0);
+    s.rc_normal       = float3(0, 1, 0);
+    s.rc_radiance     = float3(0, 0, 0);
+    s.rc_prev_pos     = primary_pos;
+    s.rc_outgoing_dir = float3(0, 1, 0);
+    s.seed_path       = seed;
+    s.path_length     = 1;
+    s.rc_length       = 2;
+    s.flags           = 0;
+    source_pdf        = 0.0f;
+
+    uint light_count = (uint)buffer_frame.restir_pt_light_count;
+    if (light_count == 0)
+        return s;
+
+    uint light_idx = min((uint)(random_float(seed) * float(light_count)), light_count - 1);
+    LightParameters light = light_parameters[light_idx];
+    if (light.intensity <= 0.0f)
+        return s;
+
+    bool is_directional = (light.flags & (1u << 0)) != 0;
+    bool is_area        = (light.flags & (1u << 6)) != 0;
+    float pick_pdf      = 1.0f / float(light_count);
+
+    if (is_directional)
+    {
+        float3 sun_dir = -light.direction;
+        float2 xi      = random_float2(seed);
+        float  cos_max = cos(SUN_CONE_HALF_ANGLE);
+        float  z       = lerp(cos_max, 1.0f, xi.x);
+        float  phi     = 2.0f * PI * xi.y;
+        float  r       = sqrt(max(1.0f - z * z, 0.0f));
+
+        float3 local = float3(cos(phi) * r, sin(phi) * r, z);
+        float3 t, b;
+        build_orthonormal_basis_fast(sun_dir, t, b);
+        float3 dir = normalize(t * local.x + b * local.y + sun_dir * z);
+
+        if (dot(dir, primary_normal) <= MIN_COS_AT_PRIMARY)
+            return s;
+
+        s.flags      |= PATH_FLAG_SKY;
+        s.rc_pos      = dir;
+        s.rc_normal   = -dir;
+        s.rc_radiance = light.color.rgb * light.intensity;
+
+        float sun_cone_pdf = 1.0f / (2.0f * PI * max(1.0f - cos_max, 1e-6f));
+        source_pdf         = pick_pdf * sun_cone_pdf;
+        return s;
+    }
+
+    if (is_area && light.area_width > 0.0f && light.area_height > 0.0f)
+    {
+        float3 light_normal = light.direction;
+        float3 light_right, light_up;
+        build_orthonormal_basis_fast(light_normal, light_right, light_up);
+
+        float2 xi = random_float2(seed);
+        float3 p  = light.position
+                  + light_right * (xi.x - 0.5f) * light.area_width
+                  + light_up    * (xi.y - 0.5f) * light.area_height;
+
+        float3 to   = p - primary_pos;
+        float  dist = length(to);
+        if (dist < 1e-3f)
+            return s;
+
+        float3 dir       = to / dist;
+        float  cos_light = dot(-dir, light_normal);
+        if (cos_light <= 0.0f || dot(dir, primary_normal) <= MIN_COS_AT_PRIMARY)
+            return s;
+
+        s.rc_pos      = p;
+        s.rc_normal   = light_normal;
+        s.rc_radiance = light.color.rgb * light.intensity;
+        s.flags      |= PATH_FLAG_HAS_RC;
+
+        float area     = light.area_width * light.area_height;
+        float area_pdf = 1.0f / area;
+        float sa_pdf   = area_pdf * dist * dist / max(cos_light, 1e-6f);
+        source_pdf     = pick_pdf * sa_pdf;
+        return s;
+    }
+
+    return s;
+}
+
 // traces a full path from the primary vertex given the first indirect direction; captures
 // x2 as the reconnection vertex candidate and accumulates the suffix radiance from x2
 // the caller samples dir via sample_brdf so the source pdf matches the primary brdf lobe
@@ -422,6 +536,7 @@ PathSample trace_path_from_primary(
     float3 primary_normal,
     float primary_roughness,
     float3 dir,
+    uint replay_seed,
     inout uint seed)
 {
     PathSample s;
@@ -430,7 +545,9 @@ PathSample trace_path_from_primary(
     s.rc_radiance     = float3(0, 0, 0);
     s.rc_prev_pos     = primary_pos;
     s.rc_outgoing_dir = float3(0, 1, 0);
-    s.seed_path       = seed;
+    // store the seed value that was used to generate xi for sample_brdf, so the random replay
+    // shift can re-derive the same xi at a destination pixel and trace a matching prefix
+    s.seed_path       = replay_seed;
     s.path_length     = 0;
     s.rc_length       = 0;
     s.flags           = 0;
@@ -464,9 +581,16 @@ PathSample trace_path_from_primary(
     s.rc_normal   = hit.geometric_normal;
     s.rc_length   = 2;
 
-    float3 suffix = accumulate_subpath_radiance(hit, -dir, RESTIR_MAX_PATH_LENGTH - 1, seed);
-    s.rc_radiance = soft_saturate_radiance(suffix, RESTIR_FIREFLY_LUMA);
-    s.path_length = 2;
+    // store the raw suffix radiance, firefly handling is deferred to the denoise temporal
+    // variance clamp and the composition stage so per-frame stochastic clamping does not
+    // flicker between frames or permanently darken bright bounces
+    float3 first_outgoing_dir;
+    float3 suffix = accumulate_subpath_radiance(hit, -dir, RESTIR_MAX_PATH_LENGTH - 1, seed, first_outgoing_dir);
+    if (any(isnan(suffix)) || any(isinf(suffix)))
+        suffix = float3(0, 0, 0);
+    s.rc_radiance     = max(suffix, 0.0f);
+    s.path_length     = 2;
+    s.rc_outgoing_dir = first_outgoing_dir;
 
     // reconnection validity: the rc vertex must be rough (stored radiance is view-dependent
     // at rc w.r.t. src's incoming, and roughness bounds the shift error) and distant enough
@@ -527,9 +651,12 @@ void ray_gen()
     // every iteration calls update_reservoir so M counts every trial (paper-form unbiased ris)
     for (uint i = 0; i < INITIAL_CANDIDATE_SAMPLES; i++)
     {
-        float2 xi = random_float2(seed);
+        // capture the seed just before consuming xi so the random replay shift can replay
+        // the same primary-direction sample at a destination pixel
+        uint   replay_seed = seed;
+        float2 xi          = random_float2(seed);
         float  source_pdf;
-        float3 dir = sample_brdf(albedo, roughness, metallic, normal_ws, view_dir, xi, source_pdf);
+        float3 dir         = sample_brdf(albedo, roughness, metallic, normal_ws, view_dir, xi, source_pdf);
 
         bool dir_valid = (source_pdf >= RESTIR_MIN_PDF) &&
                          (dot(dir, normal_ws) >= MIN_COS_AT_PRIMARY) &&
@@ -540,7 +667,7 @@ void ray_gen()
 
         if (dir_valid)
         {
-            candidate = trace_path_from_primary(pos_ws, normal_ws, roughness, dir, seed);
+            candidate = trace_path_from_primary(pos_ws, normal_ws, roughness, dir, replay_seed, seed);
             float target_pdf = target_pdf_self(candidate, pos_ws, normal_ws, view_dir, albedo, roughness, metallic);
             if (target_pdf > 0.0f)
                 weight = target_pdf / source_pdf;
@@ -549,6 +676,27 @@ void ray_gen()
         update_reservoir(reservoir, candidate, weight, random_float(seed));
     }
 
+    // additional ris stream over direct light samples (sun cone + area lights); each candidate
+    // is a single primary->light path with rc_radiance = emitted radiance and source_pdf in
+    // solid-angle measure at the primary. mixing strategies in the same reservoir keeps the
+    // estimator unbiased (each strategy independently satisfies E[weight] = integrand) and
+    // dramatically improves convergence on indirect bounces hitting concentrated emitters
+    for (uint li = 0; li < LIGHT_RIS_CANDIDATE_SAMPLES; li++)
+    {
+        float light_source_pdf;
+        PathSample light_candidate = sample_light_candidate(pos_ws, normal_ws, seed, light_source_pdf);
+
+        float light_weight = 0.0f;
+        if (light_source_pdf >= RESTIR_MIN_PDF)
+        {
+            float target_pdf = target_pdf_self(light_candidate, pos_ws, normal_ws, view_dir, albedo, roughness, metallic);
+            if (target_pdf > 0.0f)
+                light_weight = target_pdf / light_source_pdf;
+        }
+
+        update_reservoir(reservoir, light_candidate, light_weight, random_float(seed));
+    }
+
     // finalize: m_i = 1/M for the initial pass, so weight_sum becomes (1/M) * sum(p_hat/p)
     // and W = weight_sum / p_hat_y matches the paper-form output used by downstream merges
     if (reservoir.M > 0.0f)
@@ -558,10 +706,21 @@ void ray_gen()
     reservoir.target_pdf = final_target;
     reservoir.W = (final_target > 0.0f) ? (reservoir.weight_sum / final_target) : 0.0f;
 
+    // post-ris visibility test: kill the chosen sample if it is occluded so an obviously dead
+    // path does not propagate into temporal accumulation. m is preserved so the validity gate
+    // still treats this pixel as having been sampled (no spurious confidence collapse)
+    if (reservoir.W > 0.0f && !trace_shift_visibility(reservoir.sample, pos_ws, normal_ws))
+    {
+        reservoir.W           = 0.0f;
+        reservoir.weight_sum  = 0.0f;
+        reservoir.target_pdf  = 0.0f;
+    }
+
     float w_clamp = get_w_clamp_for_sample(reservoir.sample);
     reservoir.W = min(reservoir.W, w_clamp);
 
-    float sample_count_quality = saturate(reservoir.M / float(INITIAL_CANDIDATE_SAMPLES));
+    float total_candidates     = float(INITIAL_CANDIDATE_SAMPLES + LIGHT_RIS_CANDIDATE_SAMPLES);
+    float sample_count_quality = saturate(reservoir.M / max(total_candidates, 1.0f));
     reservoir.confidence       = (final_target > 0.0f) ? sample_count_quality : 0.0f;
     reservoir.age              = 0.0f;
 

data/shaders/restir_pt_denoise_spatial.hlsl

@@ -61,7 +61,9 @@ void main_cs(uint3 dispatch_id : SV_DispatchThreadID)
     float4 center_sample      = tex.Load(int3(pixel, 0));
     float3 center_color       = max(center_sample.rgb, 0.0f);
     float  center_history     = saturate(center_sample.a);
-    float  reservoir_confidence = saturate(tex_reservoir_prev4[pixel].y);
+    // confidence is the high f16 of tex4.w, see reservoir packing layout
+    uint   center_age_conf    = asuint(tex_reservoir_prev4[pixel].w);
+    float  reservoir_confidence = saturate(f16tof32(center_age_conf >> 16u));
     float  center_confidence  = saturate(center_history * 0.6f + reservoir_confidence * 0.4f);
     float  center_luma        = dot(center_color, luminance_weights);
     float  center_depth       = linearize_depth(depth);
@@ -95,7 +97,8 @@ void main_cs(uint3 dispatch_id : SV_DispatchThreadID)
         float  sample_depth   = linearize_depth(sample_depth_raw);
         float3 sample_normal  = get_normal(sample_uv);
         float4 sample_history = tex.Load(int3(sample_pixel, 0));
-        float  sample_reservoir_confidence = saturate(tex_reservoir_prev4[sample_pixel].y);
+        uint   sample_age_conf = asuint(tex_reservoir_prev4[sample_pixel].w);
+        float  sample_reservoir_confidence = saturate(f16tof32(sample_age_conf >> 16u));
         float3 sample_color   = max(sample_history.rgb, 0.0f);
         float  sample_luma    = dot(sample_color, luminance_weights);
 

data/shaders/restir_pt_denoise_temporal.hlsl

@@ -75,6 +75,83 @@ float3 clamp_history(float3 history, float3 mean, float3 sigma, float3 minimum_v
     return clamp(history, lower, upper);
 }
 
+// edge-aware 4-tap reprojection, fetches the 4 nearest history samples around prev_uv and
+// reweights them by depth + normal compatibility before bilinear blending; this stops history
+// leaking across surface edges when motion lands prev_uv between two surfaces of different depth
+float4 sample_history_edge_aware(float2 prev_uv, float3 current_normal, float current_depth, float2 history_resolution)
+{
+    float2 prev_pixel_f = prev_uv * history_resolution - 0.5f;
+    float2 base_f       = floor(prev_pixel_f);
+    float2 frac_f       = saturate(prev_pixel_f - base_f);
+
+    float4 bilinear_w = float4(
+        (1.0f - frac_f.x) * (1.0f - frac_f.y),
+        frac_f.x         * (1.0f - frac_f.y),
+        (1.0f - frac_f.x) * frac_f.y,
+        frac_f.x         * frac_f.y
+    );
+
+    int2 offsets[4] =
+    {
+        int2(0, 0), int2(1, 0), int2(0, 1), int2(1, 1)
+    };
+
+    float4 accumulated = 0.0f;
+    float  weight_sum  = 0.0f;
+    float  depth_phi   = max(current_depth * 0.05f, 1e-3f);
+
+    for (uint i = 0; i < 4; i++)
+    {
+        int2 tap_pixel = int2(base_f) + offsets[i];
+        if (tap_pixel.x < 0 || tap_pixel.x >= (int)history_resolution.x ||
+            tap_pixel.y < 0 || tap_pixel.y >= (int)history_resolution.y)
+        {
+            continue;
+        }
+
+        float2 tap_uv     = (tap_pixel + 0.5f) / history_resolution;
+        float  tap_depth_raw = tex_depth.SampleLevel(GET_SAMPLER(sampler_point_clamp), tap_uv, 0).r;
+        if (tap_depth_raw <= 0.0f)
+        {
+            continue;
+        }
+
+        float  tap_depth   = linearize_depth(tap_depth_raw);
+        float  depth_delta = abs(tap_depth - current_depth) / max(current_depth, 1e-3f);
+        if (depth_delta > 0.1f)
+        {
+            continue;
+        }
+
+        float3 tap_normal      = get_normal(tap_uv);
+        float  normal_similarity = saturate(dot(tap_normal, current_normal));
+        if (normal_similarity < 0.7f)
+        {
+            continue;
+        }
+
+        float depth_weight   = exp(-depth_delta / 0.05f);
+        float normal_weight  = pow(normal_similarity, 8.0f);
+        float w = bilinear_w[i] * depth_weight * normal_weight;
+        if (w <= 0.0f)
+        {
+            continue;
+        }
+
+        accumulated += tex2.Load(int3(tap_pixel, 0)) * w;
+        weight_sum  += w;
+    }
+
+    if (weight_sum > 0.0f)
+    {
+        return accumulated / weight_sum;
+    }
+
+    // all taps rejected, fall back to point sample at the rounded center
+    int2 fallback_pixel = clamp(int2(round(prev_pixel_f)), int2(0, 0), int2(history_resolution) - 1);
+    return tex2.Load(int3(fallback_pixel, 0));
+}
+
 bool is_history_valid(float2 current_uv, float2 prev_uv, float3 current_position, float3 current_normal, float current_depth, float2 history_resolution, out float confidence)
 {
     confidence = 0.0f;
@@ -143,7 +220,9 @@ void main_cs(uint3 dispatch_id : SV_DispatchThreadID)
     float3 current_position   = get_position(uv);
     float3 current_normal     = get_normal(uv);
     float  current_linear_z   = linearize_depth(depth);
-    float  current_confidence = saturate(tex_reservoir_prev4[pixel].y);
+    // confidence is the high f16 of tex4.w, see reservoir packing layout
+    uint   age_conf           = asuint(tex_reservoir_prev4[pixel].w);
+    float  current_confidence = saturate(f16tof32(age_conf >> 16u));
 
     float3 mean_color;
     float3 sigma_color;
@@ -164,7 +243,7 @@ void main_cs(uint3 dispatch_id : SV_DispatchThreadID)
 
     if (is_history_valid(uv, prev_uv, current_position, current_normal, current_linear_z, resolution, temporal_confidence))
     {
-        float4 history_sample = tex2.SampleLevel(GET_SAMPLER(sampler_bilinear_clamp), prev_uv, 0);
+        float4 history_sample = sample_history_edge_aware(prev_uv, current_normal, current_linear_z, float2(resolution));
         history_confidence    = saturate(history_sample.a);
         history_color         = clamp_history(max(history_sample.rgb, 0.0f), mean_color, sigma_color, minimum_color, maximum_color, current_confidence, history_confidence);