blender/intern/cycles/kernel/kernel_jitter.h
Brecht Van Lommel b20a7e01d0 Cycles: experimental correlated multi-jittered sampling pattern that can be used
instead of sobol. So far one doesn't seem to be consistently better or worse than
the other for the same number of samples but more testing is needed.

The random number generator itself is slower than sobol for most number of samples,
except 16, 64, 256, .. because they can be computed faster. This can probably be
optimized, but we can do that when/if this actually turns out to be useful.

Paper this implementation is based on:
http://graphics.pixar.com/library/MultiJitteredSampling/

Also includes some refactoring of RNG code, fixing a Sobol correlation issue with
the first BSDF and < 16 samples, skipping some unneeded RNG calls and using a
simpler unit square to unit disk function.
2013-06-07 16:06:22 +00:00

182 lines
3.6 KiB
C

/*
* Copyright 2013, Blender Foundation.
*
* This program is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License
* as published by the Free Software Foundation; either version 2
* of the License, or (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software Foundation,
* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
*/
CCL_NAMESPACE_BEGIN
/* "Correlated Multi-Jittered Sampling"
* Andrew Kensler, Pixar Technical Memo 13-01, 2013 */
/* todo: find good value, suggested 64 gives pattern on cornell box ceiling */
#define CMJ_RANDOM_OFFSET_LIMIT 4096
__device_inline bool cmj_is_pow2(int i)
{
return (i & (i - 1)) == 0;
}
__device_inline int cmj_fast_mod_pow2(int a, int b)
{
return (a & (b - 1));
}
/* a must be > 0 and b must be > 1 */
__device_inline int cmj_fast_div_pow2(int a, int b)
{
#ifdef __KERNEL_SSE2__
return a >> __builtin_ctz(b);
#else
return a/b;
#endif
}
__device_inline uint cmj_w_mask(uint w)
{
#ifdef __KERNEL_SSE2__
return ((1 << (32 - __builtin_clz(w))) - 1);
#else
w |= w >> 1;
w |= w >> 2;
w |= w >> 4;
w |= w >> 8;
w |= w >> 16;
return w;
#endif
}
__device_inline uint cmj_permute(uint i, uint l, uint p)
{
uint w = l - 1;
if((l & w) == 0) {
/* l is a power of two (fast) */
i ^= p;
i *= 0xe170893d;
i ^= p >> 16;
i ^= (i & w) >> 4;
i ^= p >> 8;
i *= 0x0929eb3f;
i ^= p >> 23;
i ^= (i & w) >> 1;
i *= 1 | p >> 27;
i *= 0x6935fa69;
i ^= (i & w) >> 11;
i *= 0x74dcb303;
i ^= (i & w) >> 2;
i *= 0x9e501cc3;
i ^= (i & w) >> 2;
i *= 0xc860a3df;
i &= w;
i ^= i >> 5;
return (i + p) & w;
}
else {
/* l is not a power of two (slow) */
w = cmj_w_mask(w);
do {
i ^= p;
i *= 0xe170893d;
i ^= p >> 16;
i ^= (i & w) >> 4;
i ^= p >> 8;
i *= 0x0929eb3f;
i ^= p >> 23;
i ^= (i & w) >> 1;
i *= 1 | p >> 27;
i *= 0x6935fa69;
i ^= (i & w) >> 11;
i *= 0x74dcb303;
i ^= (i & w) >> 2;
i *= 0x9e501cc3;
i ^= (i & w) >> 2;
i *= 0xc860a3df;
i &= w;
i ^= i >> 5;
} while (i >= l);
return (i + p) % l;
}
}
__device_inline uint cmj_hash(uint i, uint p)
{
i ^= p;
i ^= i >> 17;
i ^= i >> 10;
i *= 0xb36534e5;
i ^= i >> 12;
i ^= i >> 21;
i *= 0x93fc4795;
i ^= 0xdf6e307f;
i ^= i >> 17;
i *= 1 | p >> 18;
return i;
}
__device_inline float cmj_randfloat(uint i, uint p)
{
return cmj_hash(i, p) * (1.0f / 4294967808.0f);
}
#ifdef __CMJ__
__device_noinline float cmj_sample_1D(int s, int N, int p)
{
uint x = cmj_permute(s, N, p * 0x68bc21eb);
float jx = cmj_randfloat(s, p * 0x967a889b);
float invN = 1.0f/N;
return (x + jx)*invN;
}
__device_noinline float2 cmj_sample_2D(int s, int N, int p)
{
int m = float_to_int(sqrtf(N));
int n = (N + m - 1)/m;
float invN = 1.0f/N;
float invm = 1.0f/m;
float invn = 1.0f/n;
s = cmj_permute(s, N, p * 0x51633e2d);
int sdivm, smodm;
if(cmj_is_pow2(m)) {
sdivm = cmj_fast_div_pow2(s, m);
smodm = cmj_fast_mod_pow2(s, m);
}
else {
sdivm = float_to_int(s * invm);
smodm = s - sdivm*m;
}
uint sx = cmj_permute(smodm, m, p * 0x68bc21eb);
uint sy = cmj_permute(sdivm, n, p * 0x02e5be93);
float jx = cmj_randfloat(s, p * 0x967a889b);
float jy = cmj_randfloat(s, p * 0x368cc8b7);
return make_float2((sx + (sy + jx)*invn)*invm, (s + jy)*invN);
}
#endif
CCL_NAMESPACE_END