1//! This is a port of the original [`raymarch.hlsl`] to WGSL. It's deliberately2//! kept as close as possible so that patches to the original `raymarch.hlsl`3//! have the greatest chances of applying to this version.4//!5//! [`raymarch.hlsl`]:6//! https://gist.github.com/h3r2tic/9c8356bdaefbe80b1a22ae0aaee192db8// Copyright (c) 2023 Tomasz Stachowiak9//10// This contribution is dual licensed under EITHER OF11//12// Apache License, Version 2.0, (http://www.apache.org/licenses/LICENSE-2.0)13// MIT license (http://opensource.org/licenses/MIT)14//15// at your option.17import package::render::mesh_view_bindings::depth_prepass_texture;18import package::render::view_transformations::{19 direction_world_to_clip,20 ndc_to_uv,21 perspective_camera_near,22 position_world_to_ndc,23};25/// Allows us to sample from the depth buffer with bilinear filtering.26@if(USE_DEPTH_SAMPLERS) {27@group(2) @binding(2) var depth_linear_sampler: sampler;29/// Allows us to sample from the depth buffer with nearest-neighbor filtering.30@group(2) @binding(3) var depth_nearest_sampler: sampler;31}33/// Manual depth fetch helpers used on WebGPU where depth + filtering sampler is invalid.34@if(!USE_DEPTH_SAMPLERS)35fn depth_texel_clamped(texel: vec2<i32>) -> f32 {36 let dims = textureDimensions(depth_prepass_texture);37 let max_coord = vec2<i32>(i32(dims.x) - 1, i32(dims.y) - 1);38 let clamped = clamp(texel, vec2<i32>(0), max_coord);39 return textureLoad(depth_prepass_texture, clamped, 0);40}42@if(!USE_DEPTH_SAMPLERS)43fn depth_sample_nearest_clamped(uv: vec2<f32>, tex_size: vec2<f32>) -> f32 {44 // Match nearest sampling by snapping to the closest texel center.45 let coord = uv * tex_size - vec2(0.5);46 return depth_texel_clamped(vec2<i32>(floor(coord + vec2(0.5))));47}49@if(!USE_DEPTH_SAMPLERS)50fn depth_sample_bilinear_clamped(uv: vec2<f32>, tex_size: vec2<f32>) -> f32 {51 let coord = uv * tex_size - vec2(0.5);52 let base = vec2<i32>(floor(coord));53 let frac = coord - floor(coord);55 let d00 = depth_texel_clamped(base);56 let d10 = depth_texel_clamped(base + vec2(1, 0));57 let d01 = depth_texel_clamped(base + vec2(0, 1));58 let d11 = depth_texel_clamped(base + vec2(1, 1));60 let d0 = mix(d00, d10, frac.x);61 let d1 = mix(d01, d11, frac.x);62 return mix(d0, d1, frac.y);63}65fn depth_sample_linear(uv: vec2<f32>, tex_size: vec2<f32>) -> f32 {66 @if(USE_DEPTH_SAMPLERS)67 return textureSampleLevel(depth_prepass_texture, depth_linear_sampler, uv, 0u);68 @else69 return depth_sample_bilinear_clamped(uv, tex_size);70}72fn depth_sample_nearest(uv: vec2<f32>, tex_size: vec2<f32>) -> f32 {73 @if(USE_DEPTH_SAMPLERS)74 return textureSampleLevel(depth_prepass_texture, depth_nearest_sampler, uv, 0u);75 @else76 return depth_sample_nearest_clamped(uv, tex_size);77}79// Main code81struct HybridRootFinder {82 linear_steps: u32,83 bisection_steps: u32,84 use_secant: bool,85 linear_march_exponent: f32,87 jitter: f32,88 min_t: f32,89 max_t: f32,90}92fn hybrid_root_finder_new_with_linear_steps(v: u32) -> HybridRootFinder {93 var res: HybridRootFinder;94 res.linear_steps = v;95 res.bisection_steps = 0u;96 res.use_secant = false;97 res.linear_march_exponent = 1.0;98 res.jitter = 1.0;99 res.min_t = 0.0;100 res.max_t = 1.0;101 return res;102}104fn hybrid_root_finder_find_root(105 root_finder: ptr<function, HybridRootFinder>,106 start: vec3<f32>,107 end: vec3<f32>,108 distance_fn: ptr<function, DepthRaymarchDistanceFn>,109 hit_t: ptr<function, f32>,110 miss_t: ptr<function, f32>,111 hit_d: ptr<function, DistanceWithPenetration>,112) -> bool {113 let dir = end - start;115 var min_t = (*root_finder).min_t;116 var max_t = (*root_finder).max_t;118 var min_d = DistanceWithPenetration(0.0, false, 0.0);119 var max_d = DistanceWithPenetration(0.0, false, 0.0);121 let step_size = (max_t - min_t) / f32((*root_finder).linear_steps);123 var intersected = false;125 //126 // Ray march using linear steps128 if ((*root_finder).linear_steps > 0u) {129 let candidate_t = mix(130 min_t,131 max_t,132 pow(133 (*root_finder).jitter / f32((*root_finder).linear_steps),134 (*root_finder).linear_march_exponent135 )136 );138 let candidate = start + dir * candidate_t;139 let candidate_d = depth_raymarch_distance_fn_evaluate(distance_fn, candidate);140 intersected = candidate_d.distance < 0.0 && candidate_d.valid;142 if (intersected) {143 max_t = candidate_t;144 max_d = candidate_d;145 // The `[min_t .. max_t]` interval contains an intersection. End the linear search.146 } else {147 // No intersection yet. Carry on.148 min_t = candidate_t;149 min_d = candidate_d;151 for (var step = 1u; step < (*root_finder).linear_steps; step += 1u) {152 let candidate_t = mix(153 (*root_finder).min_t,154 (*root_finder).max_t,155 pow(156 (f32(step) + (*root_finder).jitter) / f32((*root_finder).linear_steps),157 (*root_finder).linear_march_exponent158 )159 );161 let candidate = start + dir * candidate_t;162 let candidate_d = depth_raymarch_distance_fn_evaluate(distance_fn, candidate);163 intersected = candidate_d.distance < 0.0 && candidate_d.valid;165 if (intersected) {166 max_t = candidate_t;167 max_d = candidate_d;168 // The `[min_t .. max_t]` interval contains an intersection.169 // End the linear search.170 break;171 } else {172 // No intersection yet. Carry on.173 min_t = candidate_t;174 min_d = candidate_d;175 }176 }177 }178 }180 *miss_t = min_t;181 *hit_t = min_t;183 //184 // Refine the hit using bisection186 if (intersected) {187 for (var step = 0u; step < (*root_finder).bisection_steps; step += 1u) {188 let mid_t = (min_t + max_t) * 0.5;189 let candidate = start + dir * mid_t;190 let candidate_d = depth_raymarch_distance_fn_evaluate(distance_fn, candidate);192 if (candidate_d.distance < 0.0 && candidate_d.valid) {193 // Intersection at the mid point. Refine the first half.194 max_t = mid_t;195 max_d = candidate_d;196 } else {197 // No intersection yet at the mid point. Refine the second half.198 min_t = mid_t;199 min_d = candidate_d;200 }201 }203 if ((*root_finder).use_secant) {204 // Finish with one application of the secant method205 let total_d = min_d.distance + -max_d.distance;207 let mid_t = mix(min_t, max_t, min_d.distance / total_d);208 let candidate = start + dir * mid_t;209 let candidate_d = depth_raymarch_distance_fn_evaluate(distance_fn, candidate);211 // Only accept the result of the secant method if it improves upon212 // the previous result.213 //214 // Technically root_finder should be `abs(candidate_d.distance) <215 // min(min_d.distance, -max_d.distance) * frac`, but root_finder seems216 // sufficient.217 if (abs(candidate_d.distance) < min_d.distance * 0.9 && candidate_d.valid) {218 *hit_t = mid_t;219 *hit_d = candidate_d;220 } else {221 *hit_t = max_t;222 *hit_d = max_d;223 }225 return true;226 } else {227 *hit_t = max_t;228 *hit_d = max_d;229 return true;230 }231 } else {232 // Mark the conservative miss distance.233 *hit_t = min_t;234 return false;235 }236}238struct DistanceWithPenetration {239 /// Distance to the surface of which a root we're trying to find240 distance: f32,242 /// Whether to consider this sample valid for intersection.243 /// Mostly relevant for allowing the ray marcher to travel behind surfaces,244 /// as it will mark surfaces it travels under as invalid.245 valid: bool,247 /// Conservative estimate of depth to which the ray penetrates the marched surface.248 penetration: f32,249}251struct DepthRaymarchDistanceFn {252 depth_tex_size: vec2<f32>,254 march_behind_surfaces: bool,255 depth_thickness: f32,257 use_sloppy_march: bool,258}260fn depth_raymarch_distance_fn_evaluate(261 distance_fn: ptr<function, DepthRaymarchDistanceFn>,262 ray_point_cs: vec3<f32>,263) -> DistanceWithPenetration {264 let interp_uv = ndc_to_uv(ray_point_cs.xy);266 let ray_depth = 1.0 / ray_point_cs.z;268 // We're using both point-sampled and bilinear-filtered values from the depth buffer.269 //270 // That's really stupid but works like magic. For samples taken near the ray origin,271 // the discrete nature of the depth buffer becomes a problem. It's not a land of continuous surfaces,272 // but a bunch of stacked duplo bricks.273 //274 // Technically we should be taking discrete steps in distance_fn duplo land, but then we're at the mercy275 // of arbitrary quantization of our directions -- and sometimes we'll take a step which would276 // claim that the ray is occluded -- even though the underlying smooth surface wouldn't occlude it.277 //278 // If we instead take linear taps from the depth buffer, we reconstruct the linear surface.279 // That fixes acne, but introduces false shadowing near object boundaries, as we now pretend280 // that everything is shrink-wrapped by distance_fn continuous 2.5D surface, and our depth thickness281 // heuristic ends up falling apart.282 //283 // The fix is to consider both the smooth and the discrete surfaces, and only claim occlusion284 // when the ray descends below both.285 //286 // The two approaches end up fixing each other's artifacts:287 // * The false occlusions due to duplo land are rejected because the ray stays above the smooth surface.288 // * The shrink-wrap surface is no longer continuous, so it's possible for rays to miss it.290 let linear_depth = 1.0 / depth_sample_linear(interp_uv, (*distance_fn).depth_tex_size);291 let unfiltered_depth = 1.0 / depth_sample_nearest(interp_uv, (*distance_fn).depth_tex_size);293 var max_depth: f32;294 var min_depth: f32;296 if ((*distance_fn).use_sloppy_march) {297 max_depth = unfiltered_depth;298 min_depth = unfiltered_depth;299 } else {300 max_depth = max(linear_depth, unfiltered_depth);301 min_depth = min(linear_depth, unfiltered_depth);302 }304 let bias = 0.000002;306 var res: DistanceWithPenetration;307 res.distance = max_depth * (1.0 + bias) - ray_depth;309 // distance_fn will be used at the end of the ray march to potentially discard the hit.310 res.penetration = ray_depth - min_depth;312 if ((*distance_fn).march_behind_surfaces) {313 res.valid = res.penetration < (*distance_fn).depth_thickness;314 } else {315 res.valid = true;316 }318 return res;319}321struct DepthRayMarchResult {322 /// True if the raymarch hit something.323 hit: bool,325 /// In case of a hit, the normalized distance to it.326 ///327 /// In case of a miss, the furthest the ray managed to travel, which could either be328 /// exceeding the max range, or getting behind a surface further than the depth thickness.329 ///330 /// Range: `0..=1` as a lerp factor over `ray_start_cs..=ray_end_cs`.331 hit_t: f32,333 /// UV corresponding to `hit_t`.334 hit_uv: vec2<f32>,336 /// The distance that the hit point penetrates into the hit surface.337 /// Will normally be non-zero due to limited precision of the ray march.338 ///339 /// In case of a miss: undefined.340 hit_penetration: f32,342 /// Ditto, within the range `0..DepthRayMarch::depth_thickness_linear_z`343 ///344 /// In case of a miss: undefined.345 hit_penetration_frac: f32,346}348struct DepthRayMarch {349 /// Number of steps to be taken at regular intervals to find an initial intersection.350 /// Must not be zero.351 linear_steps: u32,353 /// Exponent to be applied in the linear part of the march.354 ///355 /// A value of 1.0 will result in equidistant steps, and higher values will compress356 /// the earlier steps, and expand the later ones. This might be desirable in order357 /// to get more detail close to objects in SSR or SSGI.358 ///359 /// For optimal performance, this should be a small compile-time unsigned integer,360 /// such as 1 or 2.361 linear_march_exponent: f32,363 /// Number of steps in a bisection (binary search) to perform once the linear search364 /// has found an intersection. Helps narrow down the hit, increasing the chance of365 /// the secant method finding an accurate hit point.366 ///367 /// Useful when sampling color, e.g. SSR or SSGI, but pointless for contact shadows.368 bisection_steps: u32,370 /// Approximate the root position using the secant method -- by solving for line-line371 /// intersection between the ray approach rate and the surface gradient.372 ///373 /// Useful when sampling color, e.g. SSR or SSGI, but pointless for contact shadows.374 use_secant: bool,376 /// Jitter to apply to the first step of the linear search; 0..=1 range, mapping377 /// to the extent of a single linear step in the first phase of the search.378 /// Use 1.0 if you don't want jitter.379 jitter: f32,381 /// Clip space coordinates (w=1) of the ray.382 ray_start_cs: vec3<f32>,383 ray_end_cs: vec3<f32>,385 /// Should be used for contact shadows, but not for any color bounce, e.g. SSR.386 ///387 /// For SSR etc. this can easily create leaks, but with contact shadows it allows the rays388 /// to pass over invalid occlusions (due to thickness), and find potentially valid ones ahead.389 ///390 /// Note that this will cause the linear search to potentially miss surfaces,391 /// because when the ray overshoots and ends up penetrating a surface further than392 /// `depth_thickness_linear_z`, the ray marcher will just carry on.393 ///394 /// For this reason, this may require a lot of samples, or high depth thickness,395 /// so that `depth_thickness_linear_z >= world space ray length / linear_steps`.396 march_behind_surfaces: bool,398 /// If `true`, the ray marcher only performs nearest lookups of the depth buffer,399 /// resulting in aliasing and false occlusion when marching tiny detail.400 /// It should work fine for longer traces with fewer rays though.401 use_sloppy_march: bool,403 /// When marching the depth buffer, we only have 2.5D information, and don't know how404 /// thick surfaces are. We shall assume that the depth buffer fragments are little squares405 /// with a constant thickness defined by this parameter.406 depth_thickness_linear_z: f32,408 /// Size of the depth buffer we're marching in, in pixels.409 depth_tex_size: vec2<f32>,410}412fn depth_ray_march_new_from_depth(depth_tex_size: vec2<f32>) -> DepthRayMarch {413 var res: DepthRayMarch;414 res.jitter = 1.0;415 res.linear_steps = 4u;416 res.bisection_steps = 0u;417 res.linear_march_exponent = 1.0;418 res.depth_tex_size = depth_tex_size;419 res.depth_thickness_linear_z = 1.0;420 res.march_behind_surfaces = false;421 res.use_sloppy_march = false;422 return res;423}425fn depth_ray_march_to_cs_dir_impl(426 raymarch: ptr<function, DepthRayMarch>,427 dir_cs: vec4<f32>,428 infinite: bool,429) {430 var end_cs = vec4((*raymarch).ray_start_cs, 1.0) + dir_cs;432 // Perform perspective division, but avoid dividing by zero for rays433 // heading directly towards the eye.434 end_cs /= select(-1.0, 1.0, end_cs.w >= 0.0) * max(1e-10, abs(end_cs.w));436 // Clip ray start to the view frustum437 var delta_cs = end_cs.xyz - (*raymarch).ray_start_cs;438 let near_edge = select(vec3(-1.0, -1.0, 0.0), vec3(1.0, 1.0, 1.0), delta_cs < vec3(0.0));439 let dist_to_near_edge = (near_edge - (*raymarch).ray_start_cs) / delta_cs;440 let max_dist_to_near_edge = max(dist_to_near_edge.x, dist_to_near_edge.y);441 (*raymarch).ray_start_cs += delta_cs * max(0.0, max_dist_to_near_edge);443 // Clip ray end to the view frustum445 delta_cs = end_cs.xyz - (*raymarch).ray_start_cs;446 let far_edge = select(vec3(-1.0, -1.0, 0.0), vec3(1.0, 1.0, 1.0), delta_cs >= vec3(0.0));447 let dist_to_far_edge = (far_edge - (*raymarch).ray_start_cs) / delta_cs;448 let min_dist_to_far_edge = min(449 min(dist_to_far_edge.x, dist_to_far_edge.y),450 dist_to_far_edge.z451 );453 if (infinite) {454 delta_cs *= min_dist_to_far_edge;455 } else {456 // If unbounded, would make the ray reach the end of the frustum457 delta_cs *= min(1.0, min_dist_to_far_edge);458 }460 (*raymarch).ray_end_cs = (*raymarch).ray_start_cs + delta_cs;461}463/// March from a clip-space position (w = 1)464fn depth_ray_march_from_cs(raymarch: ptr<function, DepthRayMarch>, v: vec3<f32>) {465 (*raymarch).ray_start_cs = v;466}468/// March to a clip-space position (w = 1)469///470/// Must be called after `from_cs`, as it will clip the world-space ray to the view frustum.471fn depth_ray_march_to_cs(raymarch: ptr<function, DepthRayMarch>, end_cs: vec3<f32>) {472 let dir = vec4(end_cs - (*raymarch).ray_start_cs, 0.0) * sign(end_cs.z);473 depth_ray_march_to_cs_dir_impl(raymarch, dir, false);474}476/// March towards a clip-space direction. Infinite (ray is extended to cover the whole view frustum).477///478/// Must be called after `from_cs`, as it will clip the world-space ray to the view frustum.479fn depth_ray_march_to_cs_dir(raymarch: ptr<function, DepthRayMarch>, dir: vec4<f32>) {480 depth_ray_march_to_cs_dir_impl(raymarch, dir, true);481}483/// March to a world-space position.484///485/// Must be called after `from_cs`, as it will clip the world-space ray to the view frustum.486fn depth_ray_march_to_ws(raymarch: ptr<function, DepthRayMarch>, end: vec3<f32>) {487 depth_ray_march_to_cs(raymarch, position_world_to_ndc(end));488}490/// March towards a world-space direction. Infinite (ray is extended to cover the whole view frustum).491///492/// Must be called after `from_cs`, as it will clip the world-space ray to the view frustum.493fn depth_ray_march_to_ws_dir(raymarch: ptr<function, DepthRayMarch>, dir: vec3<f32>) {494 depth_ray_march_to_cs_dir_impl(raymarch, direction_world_to_clip(dir), true);495}497/// Perform the ray march.498fn depth_ray_march_march(raymarch: ptr<function, DepthRayMarch>) -> DepthRayMarchResult {499 var res = DepthRayMarchResult(false, 0.0, vec2(0.0), 0.0, 0.0);501 let ray_start_uv = ndc_to_uv((*raymarch).ray_start_cs.xy);502 let ray_end_uv = ndc_to_uv((*raymarch).ray_end_cs.xy);504 let ray_uv_delta = ray_end_uv - ray_start_uv;505 let ray_len_px = ray_uv_delta * (*raymarch).depth_tex_size;507 let min_px_per_step = 1u;508 let step_count = max(509 2,510 min(i32((*raymarch).linear_steps), i32(floor(length(ray_len_px) / f32(min_px_per_step))))511 );513 let linear_z_to_scaled_linear_z = 1.0 / perspective_camera_near();514 let depth_thickness = (*raymarch).depth_thickness_linear_z * linear_z_to_scaled_linear_z;516 var distance_fn: DepthRaymarchDistanceFn;517 distance_fn.depth_tex_size = (*raymarch).depth_tex_size;518 distance_fn.march_behind_surfaces = (*raymarch).march_behind_surfaces;519 distance_fn.depth_thickness = depth_thickness;520 distance_fn.use_sloppy_march = (*raymarch).use_sloppy_march;522 var hit: DistanceWithPenetration;524 var hit_t = 0.0;525 var miss_t = 0.0;526 var root_finder = hybrid_root_finder_new_with_linear_steps(u32(step_count));527 root_finder.bisection_steps = (*raymarch).bisection_steps;528 root_finder.use_secant = (*raymarch).use_secant;529 root_finder.linear_march_exponent = (*raymarch).linear_march_exponent;530 root_finder.jitter = (*raymarch).jitter;531 let intersected = hybrid_root_finder_find_root(532 &root_finder,533 (*raymarch).ray_start_cs,534 (*raymarch).ray_end_cs,535 &distance_fn,536 &hit_t,537 &miss_t,538 &hit539 );541 res.hit_t = hit_t;543 if (intersected && hit.penetration < depth_thickness && hit.distance < depth_thickness) {544 res.hit = true;545 res.hit_uv = mix(ray_start_uv, ray_end_uv, res.hit_t);546 res.hit_penetration = hit.penetration / linear_z_to_scaled_linear_z;547 res.hit_penetration_frac = hit.penetration / depth_thickness;548 return res;549 }551 res.hit_t = miss_t;552 res.hit_uv = mix(ray_start_uv, ray_end_uv, res.hit_t);554 return res;555}