073011ace31570f6fd79ca5fe07d3b6fdc4631ae
[blender-staging.git] / intern / cycles / kernel / kernel_random.h
1 /*
2  * Copyright 2011-2013 Blender Foundation
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  * http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16
17 #include "kernel/kernel_jitter.h"
18
19 CCL_NAMESPACE_BEGIN
20
21 /* Pseudo random numbers, uncomment this for debugging correlations. Only run
22  * this single threaded on a CPU for repeatable resutls. */
23 //#define __DEBUG_CORRELATION__
24
25
26 /* High Dimensional Sobol.
27  *
28  * Multidimensional sobol with generator matrices. Dimension 0 and 1 are equal
29  * to classic Van der Corput and Sobol sequences. */
30
31 #ifdef __SOBOL__
32
33 /* Skip initial numbers that are not as well distributed, especially the
34  * first sequence is just 0 everywhere, which can be problematic for e.g.
35  * path termination.
36  */
37 #define SOBOL_SKIP 64
38
39 ccl_device uint sobol_dimension(KernelGlobals *kg, int index, int dimension)
40 {
41         uint result = 0;
42         uint i = index;
43         for(uint j = 0; i; i >>= 1, j++) {
44                 if(i & 1) {
45                         result ^= kernel_tex_fetch(__sobol_directions, 32*dimension + j);
46                 }
47         }
48         return result;
49 }
50
51 #endif /* __SOBOL__ */
52
53
54 ccl_device_forceinline float path_rng_1D(KernelGlobals *kg,
55                                          RNG *rng,
56                                          int sample, int num_samples,
57                                          int dimension)
58 {
59 #ifdef __DEBUG_CORRELATION__
60         return (float)drand48();
61 #endif
62
63 #ifdef __CMJ__
64 #  ifdef __SOBOL__
65         if(kernel_data.integrator.sampling_pattern == SAMPLING_PATTERN_CMJ)
66 #  endif
67         {
68                 /* Correlated multi-jitter. */
69                 int p = *rng + dimension;
70                 return cmj_sample_1D(sample, num_samples, p);
71         }
72 #endif
73
74 #ifdef __SOBOL__
75         /* Sobol sequence value using direction vectors. */
76         uint result = sobol_dimension(kg, sample + SOBOL_SKIP, dimension);
77         float r = (float)result * (1.0f/(float)0xFFFFFFFF);
78
79         /* Cranly-Patterson rotation using rng seed */
80         float shift;
81
82         /* Hash rng with dimension to solve correlation issues.
83          * See T38710, T50116.
84          */
85         RNG tmp_rng = cmj_hash_simple(dimension, *rng);
86         shift = tmp_rng * (1.0f/(float)0xFFFFFFFF);
87
88         return r + shift - floorf(r + shift);
89 #endif
90 }
91
92 ccl_device_forceinline void path_rng_2D(KernelGlobals *kg,
93                                         RNG *rng,
94                                         int sample, int num_samples,
95                                         int dimension,
96                                         float *fx, float *fy)
97 {
98 #ifdef __DEBUG_CORRELATION__
99         *fx = (float)drand48();
100         *fy = (float)drand48();
101         return;
102 #endif
103
104 #ifdef __CMJ__
105 #  ifdef __SOBOL__
106         if(kernel_data.integrator.sampling_pattern == SAMPLING_PATTERN_CMJ)
107 #  endif
108         {
109                 /* Correlated multi-jitter. */
110                 int p = *rng + dimension;
111                 cmj_sample_2D(sample, num_samples, p, fx, fy);
112                 return;
113         }
114 #endif
115
116 #ifdef __SOBOL__
117         /* Sobol. */
118         *fx = path_rng_1D(kg, rng, sample, num_samples, dimension);
119         *fy = path_rng_1D(kg, rng, sample, num_samples, dimension + 1);
120 #endif
121 }
122
123 ccl_device_inline void path_rng_init(KernelGlobals *kg,
124                                      ccl_global uint *rng_state,
125                                      int sample, int num_samples,
126                                      RNG *rng,
127                                      int x, int y,
128                                      float *fx, float *fy)
129 {
130         /* load state */
131         *rng = *rng_state;
132         *rng ^= kernel_data.integrator.seed;
133
134 #ifdef __DEBUG_CORRELATION__
135         srand48(*rng + sample);
136 #endif
137
138         if(sample == 0) {
139                 *fx = 0.5f;
140                 *fy = 0.5f;
141         }
142         else {
143                 path_rng_2D(kg, rng, sample, num_samples, PRNG_FILTER_U, fx, fy);
144         }
145 }
146
147 /* Linear Congruential Generator */
148
149 ccl_device uint lcg_step_uint(uint *rng)
150 {
151         /* implicit mod 2^32 */
152         *rng = (1103515245*(*rng) + 12345);
153         return *rng;
154 }
155
156 ccl_device float lcg_step_float(uint *rng)
157 {
158         /* implicit mod 2^32 */
159         *rng = (1103515245*(*rng) + 12345);
160         return (float)*rng * (1.0f/(float)0xFFFFFFFF);
161 }
162
163 ccl_device uint lcg_init(uint seed)
164 {
165         uint rng = seed;
166         lcg_step_uint(&rng);
167         return rng;
168 }
169
170 /* Path Tracing Utility Functions
171  *
172  * For each random number in each step of the path we must have a unique
173  * dimension to avoid using the same sequence twice.
174  *
175  * For branches in the path we must be careful not to reuse the same number
176  * in a sequence and offset accordingly.
177  */
178
179 ccl_device_inline float path_state_rng_1D(KernelGlobals *kg,
180                                           RNG *rng,
181                                           const ccl_addr_space PathState *state,
182                                           int dimension)
183 {
184         return path_rng_1D(kg,
185                            rng,
186                            state->sample, state->num_samples,
187                            state->rng_offset + dimension);
188 }
189
190 ccl_device_inline float path_state_rng_1D_for_decision(
191         KernelGlobals *kg,
192         RNG *rng,
193         const ccl_addr_space PathState *state,
194         int dimension)
195 {
196         /* The rng_offset is not increased for transparent bounces. if we do then
197          * fully transparent objects can become subtly visible by the different
198          * sampling patterns used where the transparent object is.
199          *
200          * however for some random numbers that will determine if we next bounce
201          * is transparent we do need to increase the offset to avoid always making
202          * the same decision. */
203         const int rng_offset = state->rng_offset + state->transparent_bounce * PRNG_BOUNCE_NUM;
204         return path_rng_1D(kg,
205                            rng,
206                            state->sample, state->num_samples,
207                            rng_offset + dimension);
208 }
209
210 ccl_device_inline void path_state_rng_2D(KernelGlobals *kg,
211                                          RNG *rng,
212                                          const ccl_addr_space PathState *state,
213                                          int dimension,
214                                          float *fx, float *fy)
215 {
216         path_rng_2D(kg,
217                     rng,
218                     state->sample, state->num_samples,
219                     state->rng_offset + dimension,
220                     fx, fy);
221 }
222
223 ccl_device_inline float path_branched_rng_1D(
224         KernelGlobals *kg,
225         RNG *rng,
226         const ccl_addr_space PathState *state,
227         int branch,
228         int num_branches,
229         int dimension)
230 {
231         return path_rng_1D(kg,
232                            rng,
233                            state->sample * num_branches + branch,
234                            state->num_samples * num_branches,
235                            state->rng_offset + dimension);
236 }
237
238 ccl_device_inline float path_branched_rng_1D_for_decision(
239         KernelGlobals *kg,
240         RNG *rng,
241         const ccl_addr_space PathState *state,
242         int branch,
243         int num_branches,
244         int dimension)
245 {
246         const int rng_offset = state->rng_offset + state->transparent_bounce * PRNG_BOUNCE_NUM;
247         return path_rng_1D(kg,
248                            rng,
249                            state->sample * num_branches + branch,
250                            state->num_samples * num_branches,
251                            rng_offset + dimension);
252 }
253
254 ccl_device_inline void path_branched_rng_2D(
255         KernelGlobals *kg,
256         RNG *rng,
257         const ccl_addr_space PathState *state,
258         int branch,
259         int num_branches,
260         int dimension,
261         float *fx, float *fy)
262 {
263         path_rng_2D(kg,
264                     rng,
265                     state->sample * num_branches + branch,
266                     state->num_samples * num_branches,
267                     state->rng_offset + dimension,
268                     fx, fy);
269 }
270
271 /* Utitility functions to get light termination value,
272  * since it might not be needed in many cases.
273  */
274 ccl_device_inline float path_state_rng_light_termination(
275         KernelGlobals *kg,
276         RNG *rng,
277         const ccl_addr_space PathState *state)
278 {
279         if(kernel_data.integrator.light_inv_rr_threshold > 0.0f) {
280                 return path_state_rng_1D_for_decision(kg, rng, state, PRNG_LIGHT_TERMINATE);
281         }
282         return 0.0f;
283 }
284
285 ccl_device_inline float path_branched_rng_light_termination(
286         KernelGlobals *kg,
287         RNG *rng,
288         const ccl_addr_space PathState *state,
289         int branch,
290         int num_branches)
291 {
292         if(kernel_data.integrator.light_inv_rr_threshold > 0.0f) {
293                 return path_branched_rng_1D_for_decision(kg,
294                                                          rng,
295                                                          state,
296                                                          branch,
297                                                          num_branches,
298                                                          PRNG_LIGHT_TERMINATE);
299         }
300         return 0.0f;
301 }
302
303 ccl_device_inline void path_state_branch(ccl_addr_space PathState *state,
304                                          int branch,
305                                          int num_branches)
306 {
307         /* path is splitting into a branch, adjust so that each branch
308          * still gets a unique sample from the same sequence */
309         state->rng_offset += PRNG_BOUNCE_NUM;
310         state->sample = state->sample*num_branches + branch;
311         state->num_samples = state->num_samples*num_branches;
312 }
313
314 ccl_device_inline uint lcg_state_init(RNG *rng,
315                                       int rng_offset,
316                                       int sample,
317                                       uint scramble)
318 {
319         return lcg_init(*rng + rng_offset + sample*scramble);
320 }
321
322 ccl_device float lcg_step_float_addrspace(ccl_addr_space uint *rng)
323 {
324         /* Implicit mod 2^32 */
325         *rng = (1103515245*(*rng) + 12345);
326         return (float)*rng * (1.0f/(float)0xFFFFFFFF);
327 }
328
329 CCL_NAMESPACE_END
330