intern/cycles/kernel/closure/bsdf_ashikhmin_shirley.h

   1 /*
   2  * Copyright 2011-2014 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 __BSDF_ASHIKHMIN_SHIRLEY_H__
  18 #define __BSDF_ASHIKHMIN_SHIRLEY_H__
  19
  20 /*
  21 ASHIKHMIN SHIRLEY BSDF
  22
  23 Implementation of
  24 Michael Ashikhmin and Peter Shirley: "An Anisotropic Phong BRDF Model" (2000)
  25
  26 The Fresnel factor is missing to get a separable bsdf (intensity*color), as is
  27 the case with all other microfacet-based BSDF implementations in Cycles.
  28
  29 Other than that, the implementation directly follows the paper.
  30 */
  31
  32 CCL_NAMESPACE_BEGIN
  33
  34 ccl_device int bsdf_ashikhmin_shirley_setup(MicrofacetBsdf *bsdf)
  35 {
  36         bsdf->alpha_x = clamp(bsdf->alpha_x, 1e-4f, 1.0f);
  37         bsdf->alpha_y = bsdf->alpha_x;
  38
  39         bsdf->type = CLOSURE_BSDF_ASHIKHMIN_SHIRLEY_ID;
  40         return SD_BSDF|SD_BSDF_HAS_EVAL;
  41 }
  42
  43 ccl_device int bsdf_ashikhmin_shirley_aniso_setup(MicrofacetBsdf *bsdf)
  44 {
  45         bsdf->alpha_x = clamp(bsdf->alpha_x, 1e-4f, 1.0f);
  46         bsdf->alpha_y = clamp(bsdf->alpha_y, 1e-4f, 1.0f);
  47
  48         bsdf->type = CLOSURE_BSDF_ASHIKHMIN_SHIRLEY_ANISO_ID;
  49         return SD_BSDF|SD_BSDF_HAS_EVAL;
  50 }
  51
  52 ccl_device void bsdf_ashikhmin_shirley_blur(ShaderClosure *sc, float roughness)
  53 {
  54         MicrofacetBsdf *bsdf = (MicrofacetBsdf*)sc;
  55
  56         bsdf->alpha_x = fmaxf(roughness, bsdf->alpha_x);
  57         bsdf->alpha_y = fmaxf(roughness, bsdf->alpha_y);
  58 }
  59
  60 ccl_device_inline float bsdf_ashikhmin_shirley_roughness_to_exponent(float roughness)
  61 {
  62         return 2.0f / (roughness*roughness) - 2.0f;
  63 }
  64
  65 ccl_device_forceinline float3 bsdf_ashikhmin_shirley_eval_reflect(
  66         const ShaderClosure *sc,
  67         const float3 I,
  68         const float3 omega_in,
  69         float *pdf)
  70 {
  71         const MicrofacetBsdf *bsdf = (const MicrofacetBsdf*)sc;
  72         float3 N = bsdf->N;
  73
  74         float NdotI = dot(N, I);           /* in Cycles/OSL convention I is omega_out    */
  75         float NdotO = dot(N, omega_in);    /* and consequently we use for O omaga_in ;)  */
  76
  77         float out = 0.0f;
  78
  79         if(fmaxf(bsdf->alpha_x, bsdf->alpha_y) <= 1e-4f)
  80                 return make_float3(0.0f, 0.0f, 0.0f);
  81
  82         if(NdotI > 0.0f && NdotO > 0.0f) {
  83                 NdotI = fmaxf(NdotI, 1e-6f);
  84                 NdotO = fmaxf(NdotO, 1e-6f);
  85                 float3 H = normalize(omega_in + I);
  86                 float HdotI = fmaxf(fabsf(dot(H, I)), 1e-6f);
  87                 float HdotN = fmaxf(dot(H, N), 1e-6f);
  88
  89                 float pump = 1.0f / fmaxf(1e-6f, (HdotI*fmaxf(NdotO, NdotI))); /* pump from original paper (first derivative disc., but cancels the HdotI in the pdf nicely) */
  90                 /*float pump = 1.0f / fmaxf(1e-4f, ((NdotO + NdotI) * (NdotO*NdotI))); */ /* pump from d-brdf paper */
  91
  92                 float n_x = bsdf_ashikhmin_shirley_roughness_to_exponent(bsdf->alpha_x);
  93                 float n_y = bsdf_ashikhmin_shirley_roughness_to_exponent(bsdf->alpha_y);
  94
  95                 if(n_x == n_y) {
  96                         /* isotropic */
  97                         float e = n_x;
  98                         float lobe = powf(HdotN, e);
  99                         float norm = (n_x + 1.0f) / (8.0f * M_PI_F);
 100
 101                         out = NdotO * norm * lobe * pump;
 102                         *pdf = norm * lobe / HdotI; /* this is p_h / 4(H.I)  (conversion from 'wh measure' to 'wi measure', eq. 8 in paper) */
 103                 }
 104                 else {
 105                         /* anisotropic */
 106                         float3 X, Y;
 107                         make_orthonormals_tangent(N, bsdf->T, &X, &Y);
 108
 109                         float HdotX = dot(H, X);
 110                         float HdotY = dot(H, Y);
 111                         float lobe;
 112                         if(HdotN < 1.0f) {
 113                                 float e = (n_x * HdotX*HdotX + n_y * HdotY*HdotY) / (1.0f - HdotN*HdotN);
 114                                 lobe = powf(HdotN, e);
 115                         }
 116                         else {
 117                                 lobe = 1.0f;
 118                         }
 119                         float norm = sqrtf((n_x + 1.0f)*(n_y + 1.0f)) / (8.0f * M_PI_F);
 120
 121                         out = NdotO * norm * lobe * pump;
 122                         *pdf = norm * lobe / HdotI;
 123                 }
 124         }
 125
 126         return make_float3(out, out, out);
 127 }
 128
 129 ccl_device float3 bsdf_ashikhmin_shirley_eval_transmit(const ShaderClosure *sc, const float3 I, const float3 omega_in, float *pdf)
 130 {
 131         return make_float3(0.0f, 0.0f, 0.0f);
 132 }
 133
 134 ccl_device_inline void bsdf_ashikhmin_shirley_sample_first_quadrant(float n_x, float n_y, float randu, float randv, float *phi, float *cos_theta)
 135 {
 136         *phi = atanf(sqrtf((n_x + 1.0f) / (n_y + 1.0f)) * tanf(M_PI_2_F * randu));
 137         float cos_phi = cosf(*phi);
 138         float sin_phi = sinf(*phi);
 139         *cos_theta = powf(randv, 1.0f / (n_x * cos_phi*cos_phi + n_y * sin_phi*sin_phi + 1.0f));
 140 }
 141
 142 ccl_device int bsdf_ashikhmin_shirley_sample(const ShaderClosure *sc, float3 Ng, 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)
 143 {
 144         const MicrofacetBsdf *bsdf = (const MicrofacetBsdf*)sc;
 145         float3 N = bsdf->N;
 146         int label = LABEL_REFLECT | LABEL_GLOSSY;
 147
 148         float NdotI = dot(N, I);
 149         if(NdotI > 0.0f) {
 150
 151                 float n_x = bsdf_ashikhmin_shirley_roughness_to_exponent(bsdf->alpha_x);
 152                 float n_y = bsdf_ashikhmin_shirley_roughness_to_exponent(bsdf->alpha_y);
 153
 154                 /* get x,y basis on the surface for anisotropy */
 155                 float3 X, Y;
 156
 157                 if(n_x == n_y)
 158                         make_orthonormals(N, &X, &Y);
 159                 else
 160                         make_orthonormals_tangent(N, bsdf->T, &X, &Y);
 161
 162                 /* sample spherical coords for h in tangent space */
 163                 float phi;
 164                 float cos_theta;
 165                 if(n_x == n_y) {
 166                         /* isotropic sampling */
 167                         phi = M_2PI_F * randu;
 168                         cos_theta = powf(randv, 1.0f / (n_x + 1.0f));
 169                 }
 170                 else {
 171                         /* anisotropic sampling */
 172                         if(randu < 0.25f) {      /* first quadrant */
 173                                 float remapped_randu = 4.0f * randu;
 174                                 bsdf_ashikhmin_shirley_sample_first_quadrant(n_x, n_y, remapped_randu, randv, &phi, &cos_theta);
 175                         }
 176                         else if(randu < 0.5f) {  /* second quadrant */
 177                                 float remapped_randu = 4.0f * (.5f - randu);
 178                                 bsdf_ashikhmin_shirley_sample_first_quadrant(n_x, n_y, remapped_randu, randv, &phi, &cos_theta);
 179                                 phi = M_PI_F - phi;
 180                         }
 181                         else if(randu < 0.75f) { /* third quadrant */
 182                                 float remapped_randu = 4.0f * (randu - 0.5f);
 183                                 bsdf_ashikhmin_shirley_sample_first_quadrant(n_x, n_y, remapped_randu, randv, &phi, &cos_theta);
 184                                 phi = M_PI_F + phi;
 185                         }
 186                         else {                   /* fourth quadrant */
 187                                 float remapped_randu = 4.0f * (1.0f - randu);
 188                                 bsdf_ashikhmin_shirley_sample_first_quadrant(n_x, n_y, remapped_randu, randv, &phi, &cos_theta);
 189                                 phi = 2.0f * M_PI_F - phi;
 190                         }
 191                 }
 192
 193                 /* get half vector in tangent space */
 194                 float sin_theta = sqrtf(fmaxf(0.0f, 1.0f - cos_theta*cos_theta));
 195                 float cos_phi = cosf(phi);
 196                 float sin_phi = sinf(phi); /* no sqrt(1-cos^2) here b/c it causes artifacts */
 197                 float3 h = make_float3(
 198                         sin_theta * cos_phi,
 199                         sin_theta * sin_phi,
 200                         cos_theta
 201                         );
 202
 203                 /* half vector to world space */
 204                 float3 H = h.x*X + h.y*Y + h.z*N;
 205                 float HdotI = dot(H, I);
 206                 if(HdotI < 0.0f) H = -H;
 207
 208                 /* reflect I on H to get omega_in */
 209                 *omega_in = -I + (2.0f * HdotI) * H;
 210
 211                 if(fmaxf(bsdf->alpha_x, bsdf->alpha_y) <= 1e-4f) {
 212                         /* Some high number for MIS. */
 213                         *pdf = 1e6f;
 214                         *eval = make_float3(1e6f, 1e6f, 1e6f);
 215                         label = LABEL_REFLECT | LABEL_SINGULAR;
 216                 }
 217                 else {
 218                         /* leave the rest to eval_reflect */
 219                         *eval = bsdf_ashikhmin_shirley_eval_reflect(sc, I, *omega_in, pdf);
 220                 }
 221
 222 #ifdef __RAY_DIFFERENTIALS__
 223                 /* just do the reflection thing for now */
 224                 *domega_in_dx = (2.0f * dot(N, dIdx)) * N - dIdx;
 225                 *domega_in_dy = (2.0f * dot(N, dIdy)) * N - dIdy;
 226 #endif
 227         }
 228
 229         return label;
 230 }
 231
 232
 233 CCL_NAMESPACE_END
 234
 235 #endif /* __BSDF_ASHIKHMIN_SHIRLEY_H__ */