style cleanup: block comments
[blender.git] / intern / cycles / kernel / svm / bsdf_ward.h
1 /*
2  * Adapted from Open Shading Language with this license:
3  *
4  * Copyright (c) 2009-2010 Sony Pictures Imageworks Inc., et al.
5  * All Rights Reserved.
6  *
7  * Modifications Copyright 2011, Blender Foundation.
8  *
9  * Redistribution and use in source and binary forms, with or without
10  * modification, are permitted provided that the following conditions are
11  * met:
12  * * Redistributions of source code must retain the above copyright
13  *   notice, this list of conditions and the following disclaimer.
14  * * Redistributions in binary form must reproduce the above copyright
15  *   notice, this list of conditions and the following disclaimer in the
16  *   documentation and/or other materials provided with the distribution.
17  * * Neither the name of Sony Pictures Imageworks nor the names of its
18  *   contributors may be used to endorse or promote products derived from
19  *   this software without specific prior written permission.
20  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23  * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
24  * OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
25  * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
26  * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
27  * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
28  * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
29  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
30  * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31  */
32
33 #ifndef __BSDF_WARD_H__
34 #define __BSDF_WARD_H__
35
36 CCL_NAMESPACE_BEGIN
37
38 /* WARD */
39
40 typedef struct BsdfWardClosure {
41         //float3 m_N;
42         //float3 m_T;
43         float m_ax;
44         float m_ay;
45 } BsdfWardClosure;
46
47 __device void bsdf_ward_setup(ShaderData *sd, ShaderClosure *sc, float3 T, float ax, float ay)
48 {
49         float m_ax = clamp(ax, 1e-5f, 1.0f);
50         float m_ay = clamp(ay, 1e-5f, 1.0f);
51
52         sc->data0 = m_ax;
53         sc->data1 = m_ay;
54
55         sc->type = CLOSURE_BSDF_WARD_ID;
56         sd->flag |= SD_BSDF|SD_BSDF_HAS_EVAL|SD_BSDF_GLOSSY;
57 }
58
59 __device void bsdf_ward_blur(ShaderClosure *sc, float roughness)
60 {
61         sc->data0 = fmaxf(roughness, sc->data0);
62         sc->data1 = fmaxf(roughness, sc->data1);
63 }
64
65 __device float3 bsdf_ward_eval_reflect(const ShaderData *sd, const ShaderClosure *sc, const float3 I, const float3 omega_in, float *pdf)
66 {
67         float m_ax = sc->data0;
68         float m_ay = sc->data1;
69         float3 m_N = sd->N;
70         float3 m_T = normalize(sd->dPdu);
71
72         float cosNO = dot(m_N, I);
73         float cosNI = dot(m_N, omega_in);
74
75         if(cosNI > 0 && cosNO > 0) {
76                 // get half vector and get x,y basis on the surface for anisotropy
77                 float3 H = normalize(omega_in + I); // normalize needed for pdf
78                 float3 X, Y;
79                 make_orthonormals_tangent(m_N, m_T, &X, &Y);
80                 // eq. 4
81                 float dotx = dot(H, X) / m_ax;
82                 float doty = dot(H, Y) / m_ay;
83                 float dotn = dot(H, m_N);
84                 float exp_arg = (dotx * dotx + doty * doty) / (dotn * dotn);
85                 float denom = (4 * M_PI_F * m_ax * m_ay * sqrtf(cosNO * cosNI));
86                 float exp_val = expf(-exp_arg);
87                 float out = cosNI * exp_val / denom;
88                 float oh = dot(H, I);
89                 denom = 4 * M_PI_F * m_ax * m_ay * oh * dotn * dotn * dotn;
90                 *pdf = exp_val / denom;
91                 return make_float3 (out, out, out);
92         }
93         return make_float3 (0, 0, 0);
94 }
95
96 __device float3 bsdf_ward_eval_transmit(const ShaderData *sd, const ShaderClosure *sc, const float3 I, const float3 omega_in, float *pdf)
97 {
98         return make_float3(0.0f, 0.0f, 0.0f);
99 }
100
101 __device float bsdf_ward_albedo(const ShaderData *sd, const ShaderClosure *sc, const float3 I)
102 {
103         return 1.0f;
104 }
105
106 __device int bsdf_ward_sample(const ShaderData *sd, const ShaderClosure *sc, float randu, float randv, float3 *eval, float3 *omega_in, float3 *domega_in_dx, float3 *domega_in_dy, float *pdf)
107 {
108         float m_ax = sc->data0;
109         float m_ay = sc->data1;
110         float3 m_N = sd->N;
111         float3 m_T = normalize(sd->dPdu);
112
113         float cosNO = dot(m_N, sd->I);
114         if(cosNO > 0) {
115                 // get x,y basis on the surface for anisotropy
116                 float3 X, Y;
117                 make_orthonormals_tangent(m_N, m_T, &X, &Y);
118                 // generate random angles for the half vector
119                 // eq. 7 (taking care around discontinuities to keep
120                 //ttoutput angle in the right quadrant)
121                 // we take advantage of cos(atan(x)) == 1/sqrt(1+x^2)
122                 //tttt  and sin(atan(x)) == x/sqrt(1+x^2)
123                 float alphaRatio = m_ay / m_ax;
124                 float cosPhi, sinPhi;
125                 if(randu < 0.25f) {
126                         float val = 4 * randu;
127                         float tanPhi = alphaRatio * tanf(M_PI_2_F * val);
128                         cosPhi = 1 / sqrtf(1 + tanPhi * tanPhi);
129                         sinPhi = tanPhi * cosPhi;
130                 } else if(randu < 0.5f) {
131                         float val = 1 - 4 * (0.5f - randu);
132                         float tanPhi = alphaRatio * tanf(M_PI_2_F * val);
133                         // phi = M_PI_F - phi;
134                         cosPhi = -1 / sqrtf(1 + tanPhi * tanPhi);
135                         sinPhi = -tanPhi * cosPhi;
136                 } else if(randu < 0.75f) {
137                         float val = 4 * (randu - 0.5f);
138                         float tanPhi = alphaRatio * tanf(M_PI_2_F * val);
139                         //phi = M_PI_F + phi;
140                         cosPhi = -1 / sqrtf(1 + tanPhi * tanPhi);
141                         sinPhi = tanPhi * cosPhi;
142                 } else {
143                         float val = 1 - 4 * (1 - randu);
144                         float tanPhi = alphaRatio * tanf(M_PI_2_F * val);
145                         // phi = 2 * M_PI_F - phi;
146                         cosPhi = 1 / sqrtf(1 + tanPhi * tanPhi);
147                         sinPhi = -tanPhi * cosPhi;
148                 }
149                 // eq. 6
150                 // we take advantage of cos(atan(x)) == 1/sqrt(1+x^2)
151                 //tttt  and sin(atan(x)) == x/sqrt(1+x^2)
152                 float thetaDenom = (cosPhi * cosPhi) / (m_ax * m_ax) + (sinPhi * sinPhi) / (m_ay * m_ay);
153                 float tanTheta2 = -logf(1 - randv) / thetaDenom;
154                 float cosTheta  = 1 / sqrtf(1 + tanTheta2);
155                 float sinTheta  = cosTheta * sqrtf(tanTheta2);
156
157                 float3 h; // already normalized becaused expressed from spherical coordinates
158                 h.x = sinTheta * cosPhi;
159                 h.y = sinTheta * sinPhi;
160                 h.z = cosTheta;
161                 // compute terms that are easier in local space
162                 float dotx = h.x / m_ax;
163                 float doty = h.y / m_ay;
164                 float dotn = h.z;
165                 // transform to world space
166                 h = h.x * X + h.y * Y + h.z * m_N;
167                 // generate the final sample
168                 float oh = dot(h, sd->I);
169                 *omega_in = 2.0f * oh * h - sd->I;
170                 if(dot(sd->Ng, *omega_in) > 0) {
171                         float cosNI = dot(m_N, *omega_in);
172                         if(cosNI > 0) {
173                                 // eq. 9
174                                 float exp_arg = (dotx * dotx + doty * doty) / (dotn * dotn);
175                                 float denom = 4 * M_PI_F * m_ax * m_ay * oh * dotn * dotn * dotn;
176                                 *pdf = expf(-exp_arg) / denom;
177                                 // compiler will reuse expressions already computed
178                                 denom = (4 * M_PI_F * m_ax * m_ay * sqrtf(cosNO * cosNI));
179                                 float power = cosNI * expf(-exp_arg) / denom;
180                                 *eval = make_float3(power, power, power);
181 #ifdef __RAY_DIFFERENTIALS__
182                                 *domega_in_dx = (2 * dot(m_N, sd->dI.dx)) * m_N - sd->dI.dx;
183                                 *domega_in_dy = (2 * dot(m_N, sd->dI.dy)) * m_N - sd->dI.dy;
184                                 // Since there is some blur to this reflection, make the
185                                 // derivatives a bit bigger. In theory this varies with the
186                                 // roughness but the exact relationship is complex and
187                                 // requires more ops than are practical.
188                                 *domega_in_dx *= 10.0f;
189                                 *domega_in_dy *= 10.0f;
190 #endif
191                         }
192                 }
193         }
194         return LABEL_REFLECT|LABEL_GLOSSY;
195 }
196
197 CCL_NAMESPACE_END
198
199 #endif /* __BSDF_WARD_H__ */
200