diff --git a/Cargo.lock b/Cargo.lock
index 5b1e70c..066dc79 100644
--- a/Cargo.lock
+++ b/Cargo.lock
@@ -352,6 +352,7 @@ dependencies = [
  "openexr",
  "png_encode_mini",
  "rmath",
+ "rrand",
  "rustc-serialize",
  "scoped_threadpool",
  "sobol_burley",
@@ -565,6 +566,10 @@ dependencies = [
  "rand 0.6.5",
 ]
 
+[[package]]
+name = "rrand"
+version = "0.1.0"
+
 [[package]]
 name = "rustc-serialize"
 version = "0.3.24"
diff --git a/Cargo.toml b/Cargo.toml
index a97adbf..521ed36 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -5,6 +5,7 @@ members = [
     "sub_crates/compact",
     "sub_crates/halton",
     "sub_crates/rmath",
+    "sub_crates/rrand",
     "sub_crates/spectral_upsampling",
 ]
 
@@ -48,12 +49,15 @@ path = "sub_crates/color"
 
 [dependencies.compact]
 path = "sub_crates/compact"
-[dependencies.halton]
 
+[dependencies.halton]
 path = "sub_crates/halton"
 
 [dependencies.rmath]
 path = "sub_crates/rmath"
 
+[dependencies.rrand]
+path = "sub_crates/rrand"
+
 [dependencies.spectral_upsampling]
 path = "sub_crates/spectral_upsampling"
diff --git a/src/algorithm.rs b/src/algorithm.rs
index 36d7717..1838a9b 100644
--- a/src/algorithm.rs
+++ b/src/algorithm.rs
@@ -5,10 +5,9 @@ use std::{
     mem::MaybeUninit,
 };
 
-use crate::{
-    hash::hash_u64,
-    lerp::{lerp_slice, Lerp},
-};
+use rrand::mix_seed_u64;
+
+use crate::lerp::{lerp_slice, Lerp};
 
 /// Selects an item from a slice based on a weighting function and a
 /// number (n) between 0.0 and 1.0.  Returns the index of the selected
@@ -209,7 +208,7 @@ where
     let mut seed = n as u64;
 
     loop {
-        let i = left + (hash_u64(right as u64, seed) as usize % (right - left));
+        let i = left + (mix_seed_u64(right as u64, seed) as usize % (right - left));
 
         slc.swap(i, right - 1);
         let ii = left + {
diff --git a/src/hash.rs b/src/hash.rs
deleted file mode 100644
index 038316c..0000000
--- a/src/hash.rs
+++ /dev/null
@@ -1,45 +0,0 @@
-/// A fast seedable 32-bit hash function.
-pub fn hash_u32(mut n: u32, seed: u32) -> u32 {
-    // We rotate the bits of `seed` so it's unlikely to interact with `n`
-    // in bad ways if they're both e.g. incrementing.  The particular
-    // rotation constant used here isn't special.
-    n ^= seed.rotate_left(23);
-
-    // From https://github.com/skeeto/hash-prospector
-    n ^= n >> 16;
-    n = n.wrapping_mul(0x21f0aaad);
-    n ^= n >> 15;
-    n = n.wrapping_mul(0xd35a2d97);
-    n ^= n >> 15;
-
-    // Xor by a random number so input zero doesn't map to output zero.
-    // The particular number used here isn't special.
-    n ^ 0xe6fe3beb
-}
-
-/// A fast seedable 64-bit hash function.
-pub fn hash_u64(mut n: u64, seed: u64) -> u64 {
-    // We rotate the bits of `seed` so it's unlikely to interact with `n`
-    // in bad ways if they're both e.g. incrementing.  The particular
-    // rotation constant used here isn't special.
-    n ^= seed.rotate_left(47);
-
-    // From https://zimbry.blogspot.com/2011/09/better-bit-mixing-improving-on.html
-    n ^= n >> 30;
-    n = n.wrapping_mul(0xbf58476d1ce4e5b9);
-    n ^= n >> 27;
-    n = n.wrapping_mul(0x94d049bb133111eb);
-    n ^= n >> 31;
-
-    // Xor by a random number so input zero doesn't map to output zero.
-    // The particular number used here isn't special.
-    n ^ 0x4acc3f27cc712c9d
-}
-
-/// Returns a random float in [0, 1] based on 'n' and a seed.
-/// Generally use n for getting a bunch of different random
-/// numbers, and use seed to vary between runs.
-pub fn hash_u32_to_f32(n: u32, seed: u32) -> f32 {
-    const INV_MAX: f32 = 1.0 / std::u32::MAX as f32;
-    hash_u32(n, seed) as f32 * INV_MAX
-}
diff --git a/src/main.rs b/src/main.rs
index 45c9315..099c4ef 100644
--- a/src/main.rs
+++ b/src/main.rs
@@ -22,7 +22,6 @@ mod boundable;
 mod camera;
 mod color;
 mod fp_utils;
-mod hash;
 mod image;
 mod lerp;
 mod light;
diff --git a/src/renderer.rs b/src/renderer.rs
index 3d03950..01dec02 100644
--- a/src/renderer.rs
+++ b/src/renderer.rs
@@ -234,13 +234,16 @@ impl<'a> Renderer<'a> {
                     // the samplers themselves.
                     let si_offset =
                         owen4(morton::encode(x, y), self.seed).wrapping_mul(self.spp as u32);
-                    for si in 0..self.spp {
-                        let si = (si as u32).wrapping_add(si_offset);
+
+                    for si in 0..(self.spp as u32) {
+                        let mut sample_gen = SampleGen::new(si_offset + si, self.seed);
 
                         // Raw sample numbers.
-                        let d0 = golden_ratio_sample(si);
-                        let (d1, d2, d3, d4) = get_sample_4d(si, 0, self.seed);
-                        let (d5, _, _, _) = get_sample_4d(si, 1, self.seed);
+                        let d0 = golden_ratio_sample(si.wrapping_add(si_offset));
+                        let [d1, d2, d3, _] = sample_gen.next_dims();
+                        sample_gen.inc_padding();
+                        let [d4, d5, _, _] = sample_gen.next_dims();
+                        sample_gen.inc_padding();
 
                         // Calculate the values we need to generate a camera ray.
                         let (img_x, img_y) = {
@@ -260,7 +263,7 @@ impl<'a> Renderer<'a> {
                             .camera
                             .generate_ray(img_x, img_y, time, wavelength, lens_uv.0, lens_uv.1);
                         let path_col =
-                            trace_camera_light_path(&mut tracer, &self.scene, ray, si, self.seed);
+                            trace_camera_light_path(&mut sample_gen, &mut tracer, &self.scene, ray);
 
                         // Accummulate light path color to pixel.
                         let mut col = img_bucket.get(x, y);
@@ -322,11 +325,10 @@ impl<'a> Renderer<'a> {
 }
 
 fn trace_camera_light_path(
+    sample_gen: &mut SampleGen,
     tracer: &mut Tracer,
     scene: &Scene,
     camera_ray: Ray,
-    sample_index: u32,
-    seed: u32,
 ) -> SpectralSample {
     use crate::shading::surface_closure::SurfaceClosure;
     use crate::surface::SurfaceIntersection;
@@ -337,7 +339,6 @@ fn trace_camera_light_path(
     let mut ray_pdf = 1.0; // PDF from generating the camera/bounce ray.
     let mut acc_color = Float4::splat(0.0); // Accumulated color.
     let mut attenuation = Float4::splat(1.0); // Color attenuation along the path so far.
-    let mut sampling_seed = seed + 1;
 
     for bounce in 0..BOUNCE_COUNT {
         let isect = tracer.trace(ray);
@@ -370,8 +371,8 @@ fn trace_camera_light_path(
             //-------------------------------------------------
             // Sample light sources.
 
-            sampling_seed = sampling_seed.wrapping_add(1);
-            let (light_n, d2, d3, d4) = get_sample_4d(sample_index, 0, sampling_seed);
+            sample_gen.inc_padding();
+            let [light_n, d2, d3, d4] = sample_gen.next_dims();
 
             let light_uvw = (d2, d3, d4);
             let light_info = scene.sample_lights(
@@ -455,8 +456,8 @@ fn trace_camera_light_path(
 
             // Sample closure
             let (dir, filter, pdf) = {
-                sampling_seed = sampling_seed.wrapping_add(1);
-                let (u, v, _, _) = get_sample_4d(sample_index, 0, sampling_seed);
+                sample_gen.inc_padding();
+                let [u, v, _, _] = sample_gen.next_dims();
 
                 closure.sample(
                     idata.incoming,
@@ -509,28 +510,67 @@ fn trace_camera_light_path(
     SpectralSample::from_parts(acc_color, ray.wavelength)
 }
 
-/// Gets a sample, using LDS samples for lower dimensions,
-/// and switching to random samples at higher dimensions where
-/// LDS samples aren't available.
-#[inline(always)]
-fn get_sample_4d(i: u32, dimension_set: u32, seed: u32) -> (f32, f32, f32, f32) {
-    match dimension_set {
-        ds if ds < sobol_burley::NUM_DIMENSION_SETS_4D as u32 => {
-            // Sobol sampling.
-            //
-            let n4 = sobol_burley::sample_4d(i, ds, seed);
-            (n4[0], n4[1], n4[2], n4[3])
-        }
-        ds => {
-            // Random sampling.
-            use crate::hash::hash_u32_to_f32;
-            (
-                hash_u32_to_f32((ds * 4 + 0) ^ (i << 16), seed),
-                hash_u32_to_f32((ds * 4 + 1) ^ (i << 16), seed),
-                hash_u32_to_f32((ds * 4 + 2) ^ (i << 16), seed),
-                hash_u32_to_f32((ds * 4 + 3) ^ (i << 16), seed),
-            )
-        }
+/// Generates Owen-scrambled, padded Sobol samples.
+#[derive(Debug, Copy, Clone)]
+struct SampleGen {
+    padding_rng: rrand::Rng,
+    dimension_rng: rrand::Rng,
+    padding_seed: u32,
+    sample_i: u32,
+    sample_i_shuffled_rev: u32,
+    dimension_i: u32,
+}
+
+impl SampleGen {
+    pub fn new(sample_i: u32, seed: u32) -> Self {
+        let mut gen = Self {
+            sample_i: sample_i,
+            sample_i_shuffled_rev: 0,
+            padding_seed: 0,
+            dimension_i: 0,
+            padding_rng: rrand::Rng::new(rrand::mix_u64(seed as u64)),
+            dimension_rng: rrand::Rng::new(0),
+        };
+        gen.inc_padding();
+
+        gen
+    }
+
+    /// Increments the padding.
+    pub fn inc_padding(&mut self) {
+        use sobol_burley::parts::owen_scramble_rev;
+
+        self.padding_seed = self.padding_rng.u32();
+        self.sample_i_shuffled_rev =
+            owen_scramble_rev(self.sample_i.reverse_bits(), self.padding_seed);
+
+        self.dimension_rng = rrand::Rng::new(self.padding_rng.u64());
+        self.dimension_i = 0;
+    }
+
+    /// Gets the next four dimensions of the current sample.
+    pub fn next_dims(&mut self) -> [f32; 4] {
+        use sobol_burley::parts::{owen_scramble_int4_rev, sobol_4d_rev};
+
+        let sobol_int_rev = sobol_4d_rev(self.sample_i_shuffled_rev, self.dimension_i);
+        self.dimension_i += 1;
+
+        let rand1 = self.dimension_rng.u64();
+        let rand2 = self.dimension_rng.u64();
+
+        let sobol_owen_int = owen_scramble_int4_rev(
+            sobol_int_rev,
+            [
+                rand1 as u32,
+                (rand1 >> 32) as u32,
+                rand2 as u32,
+                (rand2 >> 32) as u32,
+            ]
+            .into(),
+        )
+        .reverse_bits();
+
+        sobol_owen_int.to_f32_norm()
     }
 }
 
diff --git a/sub_crates/rrand/Cargo.toml b/sub_crates/rrand/Cargo.toml
new file mode 100644
index 0000000..067dcf4
--- /dev/null
+++ b/sub_crates/rrand/Cargo.toml
@@ -0,0 +1,8 @@
+[package]
+name = "rrand"
+version = "0.1.0"
+edition = "2021"
+
+# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
+
+[dependencies]
diff --git a/sub_crates/rrand/src/lib.rs b/sub_crates/rrand/src/lib.rs
new file mode 100644
index 0000000..c554f2a
--- /dev/null
+++ b/sub_crates/rrand/src/lib.rs
@@ -0,0 +1,121 @@
+//! Sources of deterministic "randomness" for rendering applications.
+
+/// Convert a `u32` to a float in [0.0, 1.0).
+///
+/// Use for getting f32 values from random u32 sources.
+///
+/// Note: this is a linear mapping from [0, int_max] to [0.0, 1.0).
+#[inline(always)]
+pub fn u32_to_f32_norm(n: u32) -> f32 {
+    f32::from_bits((n >> 9) | 0x3f800000) - 1.0
+}
+
+//-------------------------------------------------------------
+
+/// A fast RNG.
+///
+#[derive(Debug, Copy, Clone)]
+pub struct Rng {
+    state: u64,
+}
+
+impl Rng {
+    /// Creates a new Rng from a seed.
+    ///
+    /// A seed of zero is perfectly fine, and does not affect the quality
+    /// of the generator.
+    #[inline]
+    pub fn new(seed: u64) -> Self {
+        Self { state: seed }
+    }
+
+    /// Gets the nth relative RNG stream from this one.
+    ///
+    /// The returned stream will be at the same point in its sequence as
+    /// this one.
+    #[inline]
+    pub fn nth_stream(&self, n: u64) -> Self {
+        Self {
+            // We just jump forward 2^40*n states.  This gives us 2^24
+            // unique streams, each of which is 2^40 numbers long.
+            state: self
+                .state
+                .wrapping_add(0xa0761d6478bd642f_u64.wrapping_mul(1 << 40).wrapping_mul(n)),
+        }
+    }
+
+    /// Returns a random u32 in [0, int_max].
+    #[inline(always)]
+    pub fn u32(&mut self) -> u32 {
+        self.u64() as u32
+    }
+
+    /// Returns a random u64 in [0, int_max].
+    #[inline(always)]
+    pub fn u64(&mut self) -> u64 {
+        // The wyrand RNG.
+        self.state = self.state.wrapping_add(0xa0761d6478bd642f);
+        let t = (self.state as u128).wrapping_mul(self.state as u128 ^ 0xe7037ed1a0b428db);
+        ((t >> 64) ^ t) as u64
+    }
+}
+
+//-------------------------------------------------------------
+
+/// A fast 32-bit mixing function.
+///
+/// Scrambles the input number to produce a different deterministic
+/// "random" number.
+#[inline(always)]
+pub fn mix_u32(mut n: u32) -> u32 {
+    // From https://github.com/skeeto/hash-prospector
+    n ^= n >> 16;
+    n = n.wrapping_mul(0x21f0aaad);
+    n ^= n >> 15;
+    n = n.wrapping_mul(0xd35a2d97);
+    n ^= n >> 15;
+
+    // Xor by a random number so input zero doesn't map to output zero.
+    // The particular number used here isn't special.
+    n ^ 0xe6fe3beb
+}
+
+/// A fast seedable 32-bit mixing function.
+///
+/// Same as `mix_u32()` but takes a seed.
+#[inline(always)]
+pub fn mix_seed_u32(n: u32, seed: u32) -> u32 {
+    // We rotate the bits of `seed` so it's unlikely to interact with `n`
+    // in bad ways if they're both e.g. incrementing.  The particular
+    // rotation constant used here isn't special.
+    mix_u32(n ^ seed.rotate_left(23))
+}
+
+/// A fast 64-bit mixing function.
+///
+/// Scrambles the input number to produce a different deterministic
+/// "random" number.
+#[inline(always)]
+pub fn mix_u64(mut n: u64) -> u64 {
+    // From https://zimbry.blogspot.com/2011/09/better-bit-mixing-improving-on.html
+    n ^= n >> 30;
+    n = n.wrapping_mul(0xbf58476d1ce4e5b9);
+    n ^= n >> 27;
+    n = n.wrapping_mul(0x94d049bb133111eb);
+    n ^= n >> 31;
+
+    // Xor by a random number so input zero doesn't map to output zero.
+    // The particular number used here isn't special.
+    n ^ 0x4acc3f27cc712c9d
+}
+
+/// A fast seedable 64-bit mixing function.
+///
+/// Same as `mix_u64()` but takes a seed.
+#[inline(always)]
+pub fn mix_seed_u64(n: u64, seed: u64) -> u64 {
+    // We rotate the bits of `seed` so it's unlikely to interact with `n`
+    // in bad ways if they're both e.g. incrementing.  The particular
+    // rotation constant used here isn't special.
+    mix_u64(n ^ seed.rotate_left(47))
+}