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