float sqr(float x) { return x*x; }
// Apply a gamma to a color
color MaxColorGamma(color c, float gamma)
{
   // Reject madness
   if (gamma <= 0.0 || gamma == 1.0)
      return c;
   // Apply a sign-safe power function
   // TODO is OSL pow sign safe ?
   return pow(c, gamma);
}
shader standard_surface(float base = .8,
                        color base_color = color(1),
                        float diffuse_roughness = 0,
                        float specular = 1,
                        color specular_color = color(1),
                        float specular_roughness = .1,
                        float specular_IOR = 1.52,
                        float specular_anisotropy = 0,
                        float specular_rotation = 0,
                        float metalness = 0,
                        float transmission = 0,
                        color transmission_color = color(1),
                        float transmission_depth = 0,
                        float transmission_scatter = 0,
                        float transmission_scatter_anisotropy = 0,
                        float transmission_dispersion = 0,
                        float transmission_extra_roughness = 0,
                        int transmit_aovs = 0
                            [[ string widget = "boolean"]],
                        float subsurface = 0,
                        color subsurface_color = color(1),
                        color subsurface_radius = color(1),
                        float subsurface_scale = 1,
                        float subsurface_anisotropy = 0,
                        int subsurface_type = 0
                            [[ string widget = "mapper",
                              string options = "diffusion:0|randomwalk:1" ]],
                        float sheen = 0,
                        color sheen_color = color(1),
                        float sheen_roughness = 0.3,
                        int thin_walled = 0
                            [[ string widget = "boolean"]],
                        normal input_normal = N,
                        vector tangent = vector(0),
                        float coat = 1,
                        color coat_color = color(1),
                        float coat_roughness = .1,
                        float coat_IOR = 1.5,
                        normal coat_normal = N,
                        float coat_affect_color = 0,
                        float coat_affect_roughness = 0,
                        float thin_film_thickness = 0,
                        float thin_film_IOR = 1.5,
                        float emission_w = 0,
                        color emission_color = color(0),
                        color opacity = color(1),
                        int caustics = 0
                            [[ string widget = "boolean"]],
                        int internal_reflections = 1
                            [[ string widget = "boolean"]],
                        int exit_to_background = 0
                            [[ string widget = "boolean"]],
                        float indirect_diffuse = 1,
                        float indirect_specular = 1,
                        output closure color result = 0)
{
   // Compute common parameters
   float metalness_clamped = clamp(metalness, 0.0, 1.0);
   float coat_clamped = clamp(coat, 0.0, 1.0);
   // Roughness
   float coat_roughness_sqr = sqr(coat_roughness);
   float rx;
   float ry;
   vector U = vector(0);
   float coat_affect = coat_clamped * coat_affect_roughness;
   float specular_roughness_sqr = 1 - (1 - sqr(specular_roughness)) * (1 - coat_affect * coat_roughness);
   if (specular_anisotropy != 0)
   {
      float aspect = sqrt(1 - clamp(specular_anisotropy, 0, 0.98));
      rx = min(specular_roughness_sqr / aspect, 1);
      ry = specular_roughness_sqr * aspect;
      // TODO: build U
   }
   else
   {
      // TODO: build U
      rx = ry = specular_roughness_sqr;
   }
   // Create closures
   closure color closures;
   // Factor to tweak color under the coat
   float coat_gamma = clamp(coat, 0.0, 1.0) * coat_affect_color + 1;
   //
   // Diffuse reflection
   //
   closures += MaxColorGamma(base * base_color, coat_gamma) * oren_nayar(input_normal, diffuse_roughness);
   //
   // Diffuse transmission
   //
   closures *= 1 - clamp(subsurface, 0, 1);
   if (thin_walled)
   {
      closures += clamp(subsurface, 0, 1) * MaxColorGamma(subsurface_color, coat_gamma) * translucent(input_normal, diffuse_roughness);
   }
   else
   {
      // Compute the final RGB mean free path
      vector subsurface_mfp = vector(subsurface_scale * subsurface_radius);
      if (subsurface_type == 0)
         closures += clamp(subsurface, 0, 1) * empirical_bssrdf(subsurface_mfp, MaxColorGamma(subsurface_color, coat_gamma));
      else
         closures += clamp(subsurface, 0, 1) * randomwalk_bssrdf(subsurface_mfp, MaxColorGamma(subsurface_color, coat_gamma), subsurface_anisotropy);
   }
   //
   // Sheen
   //
   // TODO : add Sheen closure
   //
   // Specular transmission
   //
   float transmission_IOR = (thin_walled) ? 1 : max(specular_IOR, 1.e-04);
   // TODO : handle dispersion
   // TODO : handle absorption and interior closures
   closures *= 1 - clamp(transmission, 0, 1);
   // TODO : handle specular transmission thin film
   closures += clamp(transmission, 0, 1) * transmission_color * microfacet("ggx", input_normal, U,
                                                              min(rx + sqr(transmission_extra_roughness), 1),
                                                              min(ry + sqr(transmission_extra_roughness), 1),
                                                              transmission_IOR, 1);
   //
   // Specular reflection
   //
   float Kr, Kt;
   fresnel(I, input_normal, 1 / specular_IOR, Kr, Kt);
   closures *= 1 - specular * specular_color * Kr;
   // TODO : handle handle specular reflection thin film
   closures += specular * specular_color  * microfacet("ggx", input_normal, U, rx, ry, specular_IOR, 0);
   //
   // Metal
   //
   closures *= 1 - metalness;
   // TODO : handle metal closure with complex fresnel, and thin film
   closures += metalness * base_color * microfacet("ggx", input_normal, U, rx, ry, 0, 0);
   //
   // Emission
   //
   closures += emission_w * emission_color * emission();
   //
   // Coat
   //
   // mix the coating bsdf on top of the rest based on average fresnel, tinting all
   // the BSDFs so far with the coat color
   float coat_Kr, coat_Kt;
   fresnel(I, coat_normal, 1 / coat_IOR, coat_Kr, coat_Kt);
   closures *= mix(coat, color(1), clamp(1 - Kr, 0, 1) * coat_color);
   closures += coat * microfacet("ggx", coat_normal, U, coat_roughness, coat_roughness, coat_IOR, 0);
   //
   // Opacity
   //
   closures *= opacity;
   closures += (color(1) - opacity) * transparent();
   result += closures;
}