1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
//! Provides a high quality sampling scheme based on (0, 2)-sequences
//! See sec. 7.4.3 of Physically Based Rendering

use std::{u32, f32, iter};
use rand::{Rng, StdRng};
use rand::distributions::{IndependentSample, Range};

use sampler::{Sampler, Region};

/// Low discrepancy sampler that makes use of the (0, 2) sequence to generate
/// well distributed samples
pub struct LowDiscrepancy {
    region: Region,
    /// Number of samples to take per pixel
    spp: usize,
    scramble_range: Range<u32>,
}

impl LowDiscrepancy {
    /// Create a low discrepancy sampler to sample the image in `dim.0 * dim.1` sized blocks
    pub fn new(dim: (u32, u32), mut spp: usize) -> LowDiscrepancy {
        if !spp.is_power_of_two() {
            spp = spp.next_power_of_two();
            print!("Warning: LowDiscrepancy sampler requires power of two samples per pixel, ");
            println!("rounding up to {}", spp);
        }
        LowDiscrepancy { region: Region::new((0, 0), dim), spp: spp,
                         scramble_range: Range::new(0, u32::MAX) }
    }
}

impl Sampler for LowDiscrepancy {
    fn get_samples(&mut self, samples: &mut Vec<(f32, f32)>, rng: &mut StdRng) {
        samples.clear();
        if !self.has_samples() {
            return;
        }
        if samples.len() < self.spp {
            let len = self.spp - samples.len();
            samples.extend(iter::repeat((0.0, 0.0)).take(len));
        }
        self.get_samples_2d(&mut samples[..], rng);
        for s in samples.iter_mut() {
            s.0 += self.region.current.0 as f32;
            s.1 += self.region.current.1 as f32;
        }

        self.region.current.0 += 1;
        if self.region.current.0 == self.region.end.0 {
            self.region.current.0 = self.region.start.0;
            self.region.current.1 += 1;
        }
    }
    fn get_samples_2d(&mut self, samples: &mut [(f32, f32)], rng: &mut StdRng) {
        let scramble = (self.scramble_range.ind_sample(rng),
                        self.scramble_range.ind_sample(rng));
        sample_2d(samples, scramble, 0);
        rng.shuffle(samples);
    }
    fn get_samples_1d(&mut self, samples: &mut [f32], rng: &mut StdRng) {
        let scramble = self.scramble_range.ind_sample(rng);
        sample_1d(samples, scramble, 0);
        rng.shuffle(samples);
    }
    fn max_spp(&self) -> usize { self.spp }
    fn has_samples(&self) -> bool { self.region.current.1 != self.region.end.1 }
    fn dimensions(&self) -> (u32, u32) { self.region.dim }
    fn select_block(&mut self, start: (u32, u32)) {
        self.region.select_region(start);
    }
    fn get_region(&self) -> &Region {
        &self.region
    }
}

/// Generate a 2D pattern of low discrepancy samples to fill the slice
/// sample values will be normalized between [0, 1]
pub fn sample_2d(samples: &mut [(f32, f32)], scramble: (u32, u32), offset: u32) {
    for s in samples.iter_mut().enumerate() {
        *s.1 = sample_02(s.0 as u32 + offset, scramble);
    }
}
/// Generate a 1D pattern of low discrepancy samples to fill the slice
/// sample values will be normalized between [0, 1]
pub fn sample_1d(samples: &mut [f32], scramble: u32, offset: u32) {
    for s in samples.iter_mut().enumerate() {
        *s.1 = van_der_corput(s.0 as u32 + offset, scramble);
    }
}
/// Generate a sample from a scrambled (0, 2) sequence
pub fn sample_02(n: u32, scramble: (u32, u32)) -> (f32, f32) {
    (van_der_corput(n, scramble.0), sobol(n, scramble.1))
}
/// Generate a scrambled Van der Corput sequence value
/// as described by Kollig & Keller (2002) and in PBR
/// method is specialized for base 2
pub fn van_der_corput(mut n: u32, scramble: u32) -> f32 {
	n = (n << 16) | (n >> 16);
	n = ((n & 0x00ff00ff) << 8) | ((n & 0xff00ff00) >> 8);
	n = ((n & 0x0f0f0f0f) << 4) | ((n & 0xf0f0f0f0) >> 4);
	n = ((n & 0x33333333) << 2) | ((n & 0xcccccccc) >> 2);
	n = ((n & 0x55555555) << 1) | ((n & 0xaaaaaaaa) >> 1);
	n ^= scramble;
	f32::min(((n >> 8) & 0xffffff) as f32 / ((1 << 24) as f32), 1.0 - f32::EPSILON)
}
/// Generate a scrambled Sobol' sequence value
/// as described by Kollig & Keller (2002) and in PBR
/// method is specialized for base 2
pub fn sobol(mut n: u32, mut scramble: u32) -> f32 {
    let mut i = 1 << 31;
    while n != 0 {
        if n & 0x1 != 0 {
            scramble ^= i;
        }
        n >>= 1;
        i ^= i >> 1;
    }
    f32::min(((scramble >> 8) & 0xffffff) as f32 / ((1 << 24) as f32), 1.0 - f32::EPSILON)
}