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 #ifndef __VOLUME_H__
18 #define __VOLUME_H__
19
20 CCL_NAMESPACE_BEGIN
21
22 /* VOLUME EXTINCTION */
23
volume_extinction_setup(ShaderData * sd,float3 weight)24 ccl_device void volume_extinction_setup(ShaderData *sd, float3 weight)
25 {
26 if (sd->flag & SD_EXTINCTION) {
27 sd->closure_transparent_extinction += weight;
28 }
29 else {
30 sd->flag |= SD_EXTINCTION;
31 sd->closure_transparent_extinction = weight;
32 }
33 }
34
35 /* HENYEY-GREENSTEIN CLOSURE */
36
37 typedef ccl_addr_space struct HenyeyGreensteinVolume {
38 SHADER_CLOSURE_BASE;
39
40 float g;
41 } HenyeyGreensteinVolume;
42
43 static_assert(sizeof(ShaderClosure) >= sizeof(HenyeyGreensteinVolume),
44 "HenyeyGreensteinVolume is too large!");
45
46 /* Given cosine between rays, return probability density that a photon bounces
47 * to that direction. The g parameter controls how different it is from the
48 * uniform sphere. g=0 uniform diffuse-like, g=1 close to sharp single ray. */
single_peaked_henyey_greenstein(float cos_theta,float g)49 ccl_device float single_peaked_henyey_greenstein(float cos_theta, float g)
50 {
51 return ((1.0f - g * g) / safe_powf(1.0f + g * g - 2.0f * g * cos_theta, 1.5f)) *
52 (M_1_PI_F * 0.25f);
53 };
54
volume_henyey_greenstein_setup(HenyeyGreensteinVolume * volume)55 ccl_device int volume_henyey_greenstein_setup(HenyeyGreensteinVolume *volume)
56 {
57 volume->type = CLOSURE_VOLUME_HENYEY_GREENSTEIN_ID;
58
59 /* clamp anisotropy to avoid delta function */
60 volume->g = signf(volume->g) * min(fabsf(volume->g), 1.0f - 1e-3f);
61
62 return SD_SCATTER;
63 }
64
volume_henyey_greenstein_merge(const ShaderClosure * a,const ShaderClosure * b)65 ccl_device bool volume_henyey_greenstein_merge(const ShaderClosure *a, const ShaderClosure *b)
66 {
67 const HenyeyGreensteinVolume *volume_a = (const HenyeyGreensteinVolume *)a;
68 const HenyeyGreensteinVolume *volume_b = (const HenyeyGreensteinVolume *)b;
69
70 return (volume_a->g == volume_b->g);
71 }
72
volume_henyey_greenstein_eval_phase(const ShaderClosure * sc,const float3 I,float3 omega_in,float * pdf)73 ccl_device float3 volume_henyey_greenstein_eval_phase(const ShaderClosure *sc,
74 const float3 I,
75 float3 omega_in,
76 float *pdf)
77 {
78 const HenyeyGreensteinVolume *volume = (const HenyeyGreensteinVolume *)sc;
79 float g = volume->g;
80
81 /* note that I points towards the viewer */
82 if (fabsf(g) < 1e-3f) {
83 *pdf = M_1_PI_F * 0.25f;
84 }
85 else {
86 float cos_theta = dot(-I, omega_in);
87 *pdf = single_peaked_henyey_greenstein(cos_theta, g);
88 }
89
90 return make_float3(*pdf, *pdf, *pdf);
91 }
92
93 ccl_device float3
henyey_greenstrein_sample(float3 D,float g,float randu,float randv,float * pdf)94 henyey_greenstrein_sample(float3 D, float g, float randu, float randv, float *pdf)
95 {
96 /* match pdf for small g */
97 float cos_theta;
98 bool isotropic = fabsf(g) < 1e-3f;
99
100 if (isotropic) {
101 cos_theta = (1.0f - 2.0f * randu);
102 if (pdf) {
103 *pdf = M_1_PI_F * 0.25f;
104 }
105 }
106 else {
107 float k = (1.0f - g * g) / (1.0f - g + 2.0f * g * randu);
108 cos_theta = (1.0f + g * g - k * k) / (2.0f * g);
109 if (pdf) {
110 *pdf = single_peaked_henyey_greenstein(cos_theta, g);
111 }
112 }
113
114 float sin_theta = safe_sqrtf(1.0f - cos_theta * cos_theta);
115 float phi = M_2PI_F * randv;
116 float3 dir = make_float3(sin_theta * cosf(phi), sin_theta * sinf(phi), cos_theta);
117
118 float3 T, B;
119 make_orthonormals(D, &T, &B);
120 dir = dir.x * T + dir.y * B + dir.z * D;
121
122 return dir;
123 }
124
volume_henyey_greenstein_sample(const ShaderClosure * sc,float3 I,float3 dIdx,float3 dIdy,float randu,float randv,float3 * eval,float3 * omega_in,float3 * domega_in_dx,float3 * domega_in_dy,float * pdf)125 ccl_device int volume_henyey_greenstein_sample(const ShaderClosure *sc,
126 float3 I,
127 float3 dIdx,
128 float3 dIdy,
129 float randu,
130 float randv,
131 float3 *eval,
132 float3 *omega_in,
133 float3 *domega_in_dx,
134 float3 *domega_in_dy,
135 float *pdf)
136 {
137 const HenyeyGreensteinVolume *volume = (const HenyeyGreensteinVolume *)sc;
138 float g = volume->g;
139
140 /* note that I points towards the viewer and so is used negated */
141 *omega_in = henyey_greenstrein_sample(-I, g, randu, randv, pdf);
142 *eval = make_float3(*pdf, *pdf, *pdf); /* perfect importance sampling */
143
144 #ifdef __RAY_DIFFERENTIALS__
145 /* todo: implement ray differential estimation */
146 *domega_in_dx = make_float3(0.0f, 0.0f, 0.0f);
147 *domega_in_dy = make_float3(0.0f, 0.0f, 0.0f);
148 #endif
149
150 return LABEL_VOLUME_SCATTER;
151 }
152
153 /* VOLUME CLOSURE */
154
volume_phase_eval(const ShaderData * sd,const ShaderClosure * sc,float3 omega_in,float * pdf)155 ccl_device float3 volume_phase_eval(const ShaderData *sd,
156 const ShaderClosure *sc,
157 float3 omega_in,
158 float *pdf)
159 {
160 kernel_assert(sc->type == CLOSURE_VOLUME_HENYEY_GREENSTEIN_ID);
161
162 return volume_henyey_greenstein_eval_phase(sc, sd->I, omega_in, pdf);
163 }
164
volume_phase_sample(const ShaderData * sd,const ShaderClosure * sc,float randu,float randv,float3 * eval,float3 * omega_in,differential3 * domega_in,float * pdf)165 ccl_device int volume_phase_sample(const ShaderData *sd,
166 const ShaderClosure *sc,
167 float randu,
168 float randv,
169 float3 *eval,
170 float3 *omega_in,
171 differential3 *domega_in,
172 float *pdf)
173 {
174 int label;
175
176 switch (sc->type) {
177 case CLOSURE_VOLUME_HENYEY_GREENSTEIN_ID:
178 label = volume_henyey_greenstein_sample(sc,
179 sd->I,
180 sd->dI.dx,
181 sd->dI.dy,
182 randu,
183 randv,
184 eval,
185 omega_in,
186 &domega_in->dx,
187 &domega_in->dy,
188 pdf);
189 break;
190 default:
191 *eval = make_float3(0.0f, 0.0f, 0.0f);
192 label = LABEL_NONE;
193 break;
194 }
195
196 return label;
197 }
198
199 CCL_NAMESPACE_END
200
201 #endif
202