The sampler: shuffled, Owen-scrambled Sobol (Burley 2020, as Cycles)

docs/ROADMAP.md

@@ -2318,6 +2318,13 @@ materials and the background were `@if` flags: each combination another kernel s
 merged, uniform tables, guarded and not): 0 of 190 changed; web-hello passes; Bistro's trace the same natively and
 merged.
 
+**The sampler: Owen-scrambled Sobol (2026-11-14, after it).** The kernel drew white noise from one PCG stream, its
+draws' order set by the path's branches. Now shuffled, Owen-scrambled Sobol (Burley 2020, as Cycles), each decision a
+dimension id (docs/RUNTIME.md "Random numbers"). Offline every draw is Sobol: RMSE x0.29-0.46 at 16-256 spp on the
+Cornell box, Light Room and Bistro, means within 0.01 %. Live only the camera path's first vertex, a hash elsewhere
+(as RTXPT; all-Sobol cost 2.1 ms of an 8.4 ms trace): cost and the live shots (bench-window) unchanged. The loop guard
+now skips counted loops. Next: blue-noise pixel distribution (Cycles 4.2), its own PR.
+
 ## 2. Temporal fixes
 
 Targets (2026-11-08): bench-temporal, bench-ghosting and their UE MRQ references move to Bistro night (§ 0 step 7); the

docs/RUNTIME.md

@@ -1265,6 +1265,20 @@ too is an override, from F_TINT).
 - The repository's Rust is pinned (`rust-toolchain.toml`, 1.100.0): wesl-rs 0.6 asks for 1.97.1 or newer, and the
   determinism checks' recorded hashes are of one compiler.
 
+## Random numbers
+
+`shaders/sample.wgsl`: shuffled, Owen-scrambled Sobol (Burley 2020; Cycles' Sobol-Burley), four dimensions padded by
+dimension ids. Each decision at a path vertex has an id (`rng_dim`: path, vertex, slot `D_*`; `rng_sub`: a loop's
+iteration or a light), which seeds the scramble: one id's samples are stratified over the sample index (the sample
+offline, the frame live; 16 bits), ids are uncorrelated.
+- Offline: Sobol for every draw. Live: only the camera path's first vertex; loops, later vertices and the training
+  paths use a hash, as NVIDIA's RTXPT (all-Sobol cost 2.1 ms of Bistro's 8.4 ms live trace; this costs nothing
+  measurable).
+- Offline RMSE against the other sampler at 4096 spp (320 x 180, no denoiser): x0.29-0.46 at 16-256 spp (Cornell
+  box, Light Room, Bistro), x1 at 1 spp; means within 0.01 %. `sample.rs` checks the (0, m, 2)-net property and the
+  shader against itself on the GPU.
+- The loop guard leaves counted loops (`for (var i = N; i < M; i++)`, literals) alone, as Tint does.
+
 ## What the kernel trusts
 
 The kernel follows indices it reads from buffers and loops over counts it is given; it checks none of them (a

renderer/shaders/sample.wgsl

@@ -0,0 +1,176 @@
+// The kernel's random numbers: shuffled, Owen-scrambled Sobol (Burley 2020; Cycles' Sobol-Burley). Four dimensions,
+// more by padding: each decision has a dimension id that seeds the scramble. renderer/src/sample.rs: the same on the CPU.
+
+// Sobol dimensions 1-3 (Joe & Kuo, bit-reversed) of a 16-bit bit-reversed index, unrolled. Dimension 0: reverseBits.
+fn sobol_d1(x: u32) -> u32 {
+    var r = 0u;
+    r ^= select(0u, 0x00000001u, (x & 0x80000000u) != 0u);
+    r ^= select(0u, 0x00000003u, (x & 0x40000000u) != 0u);
+    r ^= select(0u, 0x00000005u, (x & 0x20000000u) != 0u);
+    r ^= select(0u, 0x0000000fu, (x & 0x10000000u) != 0u);
+    r ^= select(0u, 0x00000011u, (x & 0x08000000u) != 0u);
+    r ^= select(0u, 0x00000033u, (x & 0x04000000u) != 0u);
+    r ^= select(0u, 0x00000055u, (x & 0x02000000u) != 0u);
+    r ^= select(0u, 0x000000ffu, (x & 0x01000000u) != 0u);
+    r ^= select(0u, 0x00000101u, (x & 0x00800000u) != 0u);
+    r ^= select(0u, 0x00000303u, (x & 0x00400000u) != 0u);
+    r ^= select(0u, 0x00000505u, (x & 0x00200000u) != 0u);
+    r ^= select(0u, 0x00000f0fu, (x & 0x00100000u) != 0u);
+    r ^= select(0u, 0x00001111u, (x & 0x00080000u) != 0u);
+    r ^= select(0u, 0x00003333u, (x & 0x00040000u) != 0u);
+    r ^= select(0u, 0x00005555u, (x & 0x00020000u) != 0u);
+    r ^= select(0u, 0x0000ffffu, (x & 0x00010000u) != 0u);
+    return r;
+}
+fn sobol_d2(x: u32) -> u32 {
+    var r = 0u;
+    r ^= select(0u, 0x00000001u, (x & 0x80000000u) != 0u);
+    r ^= select(0u, 0x00000003u, (x & 0x40000000u) != 0u);
+    r ^= select(0u, 0x00000006u, (x & 0x20000000u) != 0u);
+    r ^= select(0u, 0x00000009u, (x & 0x10000000u) != 0u);
+    r ^= select(0u, 0x00000017u, (x & 0x08000000u) != 0u);
+    r ^= select(0u, 0x0000003au, (x & 0x04000000u) != 0u);
+    r ^= select(0u, 0x00000071u, (x & 0x02000000u) != 0u);
+    r ^= select(0u, 0x000000a3u, (x & 0x01000000u) != 0u);
+    r ^= select(0u, 0x00000116u, (x & 0x00800000u) != 0u);
+    r ^= select(0u, 0x00000339u, (x & 0x00400000u) != 0u);
+    r ^= select(0u, 0x00000677u, (x & 0x00200000u) != 0u);
+    r ^= select(0u, 0x000009aau, (x & 0x00100000u) != 0u);
+    r ^= select(0u, 0x00001601u, (x & 0x00080000u) != 0u);
+    r ^= select(0u, 0x00003903u, (x & 0x00040000u) != 0u);
+    r ^= select(0u, 0x00007706u, (x & 0x00020000u) != 0u);
+    r ^= select(0u, 0x0000aa09u, (x & 0x00010000u) != 0u);
+    return r;
+}
+fn sobol_d3(x: u32) -> u32 {
+    var r = 0u;
+    r ^= select(0u, 0x00000001u, (x & 0x80000000u) != 0u);
+    r ^= select(0u, 0x00000003u, (x & 0x40000000u) != 0u);
+    r ^= select(0u, 0x00000004u, (x & 0x20000000u) != 0u);
+    r ^= select(0u, 0x0000000au, (x & 0x10000000u) != 0u);
+    r ^= select(0u, 0x0000001fu, (x & 0x08000000u) != 0u);
+    r ^= select(0u, 0x0000002eu, (x & 0x04000000u) != 0u);
+    r ^= select(0u, 0x00000045u, (x & 0x02000000u) != 0u);
+    r ^= select(0u, 0x000000c9u, (x & 0x01000000u) != 0u);
+    r ^= select(0u, 0x0000011bu, (x & 0x00800000u) != 0u);
+    r ^= select(0u, 0x000002a4u, (x & 0x00400000u) != 0u);
+    r ^= select(0u, 0x0000079au, (x & 0x00200000u) != 0u);
+    r ^= select(0u, 0x00000b67u, (x & 0x00100000u) != 0u);
+    r ^= select(0u, 0x0000101eu, (x & 0x00080000u) != 0u);
+    r ^= select(0u, 0x0000302du, (x & 0x00040000u) != 0u);
+    r ^= select(0u, 0x00004041u, (x & 0x00020000u) != 0u);
+    r ^= select(0u, 0x0000a0c3u, (x & 0x00010000u) != 0u);
+    return r;
+}
+
+// Owen scrambling of a bit-reversed integer (Vegdahl 2021's hash, a better Laine-Karras permutation).
+fn owen_rev(n_in: u32, seed: u32) -> u32 {
+    var n = n_in;
+    n ^= n * 0x3d20adeau;
+    n += seed;
+    n *= (seed >> 16u) | 1u;
+    n ^= n * 0x05526c56u;
+    n ^= n * 0x53a22864u;
+    return n;
+}
+
+fn hash_hp(x: u32) -> u32 {
+    var i = x;
+    i ^= i >> 16u;
+    i *= 0x21f0aaadu;
+    i ^= i >> 15u;
+    i *= 0xd35a2d97u;
+    i ^= i >> 15u;
+    return i ^ 0xe6fe3bebu;
+}
+
+// One dimension (0-3) of the Sobol point at a bit-reversed index, Owen-scrambled.
+fn sobol_owen(rev_index: u32, dim: u32, seed: u32) -> f32 {
+    var r = 0u;
+    if (dim == 0u) {
+        r = reverseBits(rev_index);
+    } else if (dim == 1u) {
+        r = sobol_d1(rev_index);
+    } else if (dim == 2u) {
+        r = sobol_d2(rev_index);
+    } else {
+        r = sobol_d3(rev_index);
+    }
+    r = reverseBits(owen_rev(r, seed));
+    return f32(r) * (1.0 / 4294967808.0);
+}
+
+// A path's sample: its index in the pixel's sequence, the pixel's seed, which path of the invocation it is.
+var<private> RNG_INDEX: u32;
+var<private> RNG_SEED: u32;
+var<private> RNG_PATH: u32;
+
+// The id of decision `slot` at the path's current vertex (`cur_depth`); `rng_sub` the k-th part of decision `d` (a
+// loop's iteration, a light).
+fn rng_dim(slot: u32) -> u32 { return (RNG_PATH << 28u) ^ (cur_depth << 20u) ^ slot; }
+fn rng_sub(d: u32, k: u32) -> u32 { return hash_hp(d ^ hash_hp(k + 0x9e3779b9u)) | 0x80000000u; }
+
+// A path vertex's decisions (`rng_dim`).
+const D_CAMERA: u32 = 0u;    // pixel jitter, lens (4D)
+const D_BSDF: u32 = 1u;      // lobe, specular or not, direction (4D)
+const D_RR: u32 = 2u;        // Russian roulette
+const D_PICK: u32 = 3u;      // light_pick's kind of light
+const D_SUN: u32 = 4u;       // a directional light's cone (2D a light)
+const D_SUN_VIS: u32 = 5u;   // its shadow ray
+const D_LOCAL: u32 = 6u;     // local lights (`local_lights`)
+const D_EMIT: u32 = 7u;      // emitters (`sample_emitter_ris`, `sample_emitter`)
+const D_EMIT_VIS: u32 = 8u;
+const D_ENV: u32 = 9u;       // the environment map (4D)
+const D_ENV_VIS: u32 = 10u;
+const D_FOG: u32 = 11u;      // height fog's in-scattering (`fog_scatter`)
+const D_MEDIA: u32 = 12u;    // bounded media (`media_track`, `media_march`)
+const D_PHASE: u32 = 13u;    // a medium's scattered direction (`sample_hg`)
+const D_VOL_LIGHT: u32 = 14u; // lights at a point in a medium (`phase_direct`, `point_direct`)
+const D_VOXEL: u32 = 15u;    // the radiance cache's volume points (`point_sphere`)
+const D_CACHE: u32 = 16u;    // the cache's fixed-point dither
+
+// Live: Sobol only for the camera path's first vertex (not a loop's parts: `rng_sub` sets the top bit), a hash
+// elsewhere (as NVIDIA's RTXPT; all-Sobol cost 2 ms of Bistro's trace). Offline: Sobol everywhere (as Cycles).
+fn rng_white(d: u32) -> bool { return LIVE && (cur_depth > 0u || RNG_PATH != 0u || d >= 0x80000000u); }
+fn white(d: u32, c: u32) -> f32 { return f32(hash_hp(rng_seed(d) ^ (RNG_INDEX * 0x9e3779b9u + c * 0x85ebca6bu)) >> 8u) * (1.0 / 16777216.0); }
+
+// A hash, always: for loops of no fixed length (delta tracking's steps), where a new Sobol dimension a step is QMC's
+// worst case (as pbrt-v4, which draws its majorant tracking from a hashed RNG).
+fn rng_hash(d: u32, c: u32) -> f32 { return white(d, c); }
+fn rng_hash3(d: u32) -> vec3<f32> { return vec3<f32>(white(d, 0u), white(d, 1u), white(d, 2u)); }
+
+// Dimension id `d`'s scramble seed; the shuffled index, 16 bits (65,536 samples, then it repeats; Cycles' index mask).
+fn rng_seed(d: u32) -> u32 { return RNG_SEED ^ hash_hp(d); }
+fn rng_index(salt: u32) -> u32 { return owen_rev(reverseBits(RNG_INDEX), salt) & 0xffff0000u; }
+
+fn rng1(d: u32) -> f32 {
+    if (rng_white(d)) { return white(d, 0u); }
+    let s = rng_seed(d);
+    return sobol_owen(rng_index(s ^ 0xbff95bfeu), 0u, s ^ 0x635c77bdu);
+}
+fn rng2(d: u32) -> vec2<f32> {
+    if (rng_white(d)) { return vec2<f32>(white(d, 0u), white(d, 1u)); }
+    let s = rng_seed(d);
+    let i = rng_index(s ^ 0xf8ade99au);
+    return vec2<f32>(sobol_owen(i, 0u, s ^ 0xe0aaaf76u), sobol_owen(i, 1u, s ^ 0x94964d4eu));
+}
+fn rng3(d: u32) -> vec3<f32> {
+    if (rng_white(d)) { return vec3<f32>(white(d, 0u), white(d, 1u), white(d, 2u)); }
+    let s = rng_seed(d);
+    let i = rng_index(s ^ 0xcaa726acu);
+    return vec3<f32>(sobol_owen(i, 0u, s ^ 0x9e78e391u), sobol_owen(i, 1u, s ^ 0x67c33241u), sobol_owen(i, 2u, s ^ 0x78c395c5u));
+}
+fn rng4(d: u32) -> vec4<f32> {
+    if (rng_white(d)) { return vec4<f32>(white(d, 0u), white(d, 1u), white(d, 2u), white(d, 3u)); }
+    let s = rng_seed(d);
+    let i = rng_index(s ^ 0xc2c1a055u);
+    return vec4<f32>(sobol_owen(i, 0u, s ^ 0x39468210u), sobol_owen(i, 1u, s ^ 0xe9d8a845u), sobol_owen(i, 2u, s ^ 0x5f32b482u), sobol_owen(i, 3u, s ^ 0x1524cc56u));
+}
+
+// Starts sample `index` of a pixel (the sample offline, the frame live) on path `path`; `seed`: the render's
+// (THREERT_SEED offline) and the stream's.
+fn rng_start(pix: vec2<u32>, index: u32, seed: u32, path: u32) {
+    RNG_SEED = hash_hp(pix.x ^ hash_hp(pix.y ^ hash_hp(seed)));
+    RNG_INDEX = index;
+    RNG_PATH = path;
+}

renderer/shaders/translucent.wgsl

@@ -84,15 +84,16 @@ fn translucent(@builtin(global_invocation_id) gid: vec3<u32>) {
     let h = select(TPS.size.y, TPS.view.y, low);
     if (gid.x >= w || gid.y >= h) { return; }
     let pix = gid.y * w + gid.x;
-    rng_state = pcg(pix ^ pcg(P.size.z * 9781u + P.size.w * 6271u + 0x5bd1e995u));
+    rng_start(gid.xy, P.size.z, P.size.w ^ 0x5bd1e995u, 4u);
     cur_depth = 0u;
+    let cam = rng2(rng_dim(D_CAMERA));
     mirror_path = false;
     // Live: the pixel centre (stable streaks; the image has passed the temporal filter). Offline: a random point in it.
     var jx = 0.5;
     var jy = 0.5;
     if (!live) {
-        jx = rnd();
-        jy = rnd();
+        jx = cam.x;
+        jy = cam.y;
     }
     let ndc = vec2<f32>((f32(gid.x) + jx) / f32(w) * 2.0 - 1.0, 1.0 - (f32(gid.y) + jy) / f32(h) * 2.0);
     let d = normalize(P.cam_fwd.xyz + P.cam_right.xyz * (ndc.x * P.cam_pos.w * P.cam_fwd.w) + P.cam_up.xyz * (ndc.y * P.cam_pos.w));
@@ -140,7 +141,7 @@ fn translucent(@builtin(global_invocation_id) gid: vec3<u32>) {
     // Front to back.
     var col = vec3<f32>(0.0);
     var tr = 1.0;
-    var train = live && rnd() < TPS.misc.x;
+    var train = live && rng1(rng_dim(D_CACHE)) < TPS.misc.x;
     let foot = 2.0 * P.cam_pos.w / f32(h);
     for (var k = 0u; k < n; k++) {
         let ti = TPS.size.z + lp[k];
@@ -203,7 +204,7 @@ fn translucent_up(@builtin(global_invocation_id) gid: vec3<u32>) {
 fn volume_train(@builtin(global_invocation_id) gid: vec3<u32>) {
     let k = gid.x;
     if (k >= min(atomicLoad(&vcount[0]), TPS.view.z)) { return; }
-    rng_state = pcg(k ^ pcg(P.size.z * 7919u + P.size.w * 104729u + 0x68e31da4u));
+    rng_start(vec2<u32>(k, 0u), P.size.z, P.size.w ^ 0x68e31da4u, 5u);
     mirror_path = false;
     let x = vpts[k].xyz;
     let e = point_radiance_train(x, P.misc.x);

renderer/src/gpu.rs

@@ -698,8 +698,9 @@ pub fn instance_mask(cast_shadow: bool, hidden: bool) -> u8 {
 }
 
 /// The kernel's files, one module in this order (pathtrace.wgsl's head says what each holds and names the `@if` flags).
-const KERNEL_FILES: [&str; 4] = [
+const KERNEL_FILES: [&str; 5] = [
     include_str!("../shaders/pathtrace.wgsl"),
+    include_str!("../shaders/sample.wgsl"),
     include_str!("../shaders/trace_units.wgsl"),
     include_str!("../shaders/trace_hw.wgsl"),
     include_str!("../shaders/translucent.wgsl"),
@@ -821,6 +822,34 @@ impl wesl::resolver::Resolver for KernelResolver<'_> {
 /// scene's functions' too, begins with `if (loop_guard()) { break; }` (pathtrace.wgsl declares `loop_guard`, under
 /// `@if(guard)`). A pass over the variant's syntax tree, as a compiler instruments code: the tree knows what a loop
 /// and its body are, so nothing the text says (comments, strings, a loop without braces) can mislead it.
+/// A `for` loop that ends by its own count: `for (var i = N; i < M; i++)` (or `<=`), N and M integer literals, whose
+/// body does not assign `i`. The guard leaves it out, as Tint's PreventInfiniteLoops leaves loops it can prove finite.
+fn counted(f: &wesl::syntax::ForStatement) -> bool {
+    let (Some(init), Some(cond), Some(update)) = (&f.initializer, &f.condition, &f.update) else { return false };
+    let lit = |t: &str| {
+        let t = t.trim().trim_end_matches(['u', 'i']);
+        !t.is_empty() && t.bytes().all(|c| c.is_ascii_digit())
+    };
+    let init = init.to_string();
+    let Some((name, start)) = init.trim().trim_end_matches(';').strip_prefix("var ").and_then(|r| r.split_once('=')) else { return false };
+    let name = name.split(':').next().unwrap_or("").trim();
+    let cond = cond.to_string();
+    let bound = cond.trim().trim_start_matches('(').trim_end_matches(')').strip_prefix(name).map(str::trim_start).and_then(|r| r.strip_prefix("<=").or_else(|| r.strip_prefix('<')));
+    if name.is_empty() || !lit(start) || !bound.is_some_and(lit) || update.to_string().trim().trim_end_matches(';') != format!("{name}++") {
+        return false;
+    }
+    // Not assigned in the body: `name` as an identifier followed by `=` (not `==`), `op=`, `++` or `--`.
+    let b = f.body.to_string();
+    let bytes = b.as_bytes();
+    let ident = |c: u8| c.is_ascii_alphanumeric() || c == b'_';
+    !b.match_indices(name).any(|(k, _)| {
+        let before = k.checked_sub(1).map(|p| bytes[p]);
+        let after = b[k + name.len()..].trim_start();
+        before.is_none_or(|c| !ident(c) && c != b'.') && !after.starts_with(|c: char| c.is_ascii_alphanumeric() || c == '_')
+            && (after.starts_with("++") || after.starts_with("--") || (after.starts_with('=') && !after.starts_with("==")) || ["+=", "-=", "*=", "/=", "%=", "&=", "|=", "^=", "<<=", ">>="].iter().any(|op| after.starts_with(op)))
+    })
+}
+
 pub fn guard_loops(tu: &mut wesl::syntax::TranslationUnit) {
     use wesl::syntax::{CompoundStatement, GlobalDeclaration, Statement, StatementNode};
     fn body(b: &mut CompoundStatement, guard: &StatementNode) {
@@ -851,6 +880,7 @@ pub fn guard_loops(tu: &mut wesl::syntax::TranslationUnit) {
                     body(&mut c.body, guard);
                 }
             }
+            Statement::For(f) if counted(f) => body(&mut f.body, guard),
             Statement::For(f) => opened(&mut f.body, guard),
             Statement::While(w) => opened(&mut w.body, guard),
             _ => {}

renderer/src/lib.rs

@@ -23,6 +23,7 @@ pub mod scene;
 #[cfg(not(target_arch = "wasm32"))]
 pub mod script;
 pub mod skin;
+pub mod sample;
 pub mod sky;
 pub mod tables;
 pub mod units;